{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "858c5f0b",
   "metadata": {},
   "source": [
    "# LU Decompostion\n",
    "**강좌**: *수치해석*"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "57a6f719",
   "metadata": {},
   "source": [
    "## LU 분해 개요\n",
    "- Gauss Elimination은 선형 방정식 $Ax=b$ 를 Upper triangular matrix $U$를 이용한 $Ux=c$ 로 바꿔서 푸는 방법이다. \n",
    "- 이론적으로 $A=LU$ 로 분해해서 푸는 것과 같다.\n",
    "\n",
    "$$\n",
    "    Ax=LUx=b\\\\\n",
    "    Ux=c, Lc=b\n",
    "$$\n",
    "\n",
    "- $A=LU$ 로 분해해놓으면 $b$ 벡터가 달라져도 Substitution 과정만으로 빠르게 계산할 수 있다."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "916ae787",
   "metadata": {},
   "source": [
    "## 이론\n",
    "- Gauss Elimination 과정을 Matrix Operation으로 생각하면 다음과 같다.\n",
    "    * 아래 예제를 생각하자.\n",
    "        \n",
    "        $$\n",
    "        A = \\begin{bmatrix}\n",
    "        2 & 1 & 1 \\\\\n",
    "        4 & 6 & 0 \\\\\n",
    "        -2 & 7 & 2\n",
    "        \\end{bmatrix}\n",
    "        $$\n",
    "        \n",
    "    * $a_{21}$ 을 제거하기 위한 Row Operation ($(2) - 2\\times(1)$) 연산 Matrix\n",
    "    \n",
    "        $$\n",
    "        E_{21}=\\left[\\begin{array}{ccc}\n",
    "        1 & 0 & 0\\\\\n",
    "        -2 & 1 & 0\\\\\n",
    "        0 & 0 & 1\n",
    "        \\end{array}\\right]$$\n",
    "        \n",
    "    * $a_{31}$ 을 제거하기 위한 Row Operation ($(3) + 1\\times(1)$) 연산 Matrix\n",
    "    \n",
    "        $$E_{31}=\\left[\\begin{array}{ccc}\n",
    "        1 & 0 & 0\\\\\n",
    "        0 & 1 & 0\\\\\n",
    "        l & 0 & 1\n",
    "        \\end{array}\\right]$$       \n",
    "        \n",
    "    * $a_{32}$ 을 제거하기 위한 Row Operation ($(2) + 1\\times(1)$) 연산 Matrix\n",
    "    \n",
    "        $$E_{32}=\\left[\\begin{array}{ccc}\n",
    "        1 & 0 & 0\\\\\n",
    "        0 & 1 & 0\\\\\n",
    "        0 & l & 1\n",
    "        \\end{array}\\right]$$       "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "b51cb2a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "\n",
    "def ematrix(n,i,j,d):\n",
    "    \"\"\"\n",
    "    Elimination matrix\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    n : integer\n",
    "        size of matrix\n",
    "    i : integer\n",
    "        row of pivot\n",
    "    j : intger\n",
    "        row of elimination\n",
    "    d : float\n",
    "        elimination coefficeint\n",
    "        \n",
    "    Returns\n",
    "    -------\n",
    "    mat : array\n",
    "        elimination matrix\n",
    "    \"\"\"\n",
    "    # Make Identity matrix (size n)\n",
    "    mat = np.eye(n)        \n",
    "    \n",
    "    # Assign elimination coefficient\n",
    "    mat[j,i] = d\n",
    "    \n",
    "    return mat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "d0e09b1e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Problem\n",
    "A = np.array([[2,1,1],[4,-6,0],[-2,7,2]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "495dbf38",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2.  1.  1.]\n",
      " [ 0. -8. -2.]\n",
      " [-2.  7.  2.]]\n"
     ]
    }
   ],
   "source": [
    "# One step (E21)\n",
    "E21 = ematrix(3, 0, 1, -2)\n",
    "print(E21 @ A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "dc2d9dd7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2.  1.  1.]\n",
      " [ 0. -8. -2.]\n",
      " [ 0.  0.  1.]]\n"
     ]
    }
   ],
   "source": [
    "# For all steps\n",
    "E31 = ematrix(3, 0, 2, 1)\n",
    "E32 = ematrix(3, 1, 2, 1)\n",
    "U = E32 @ E31 @ E21 @ A \n",
    "print(U)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "4ad76630",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.,  0.,  0.],\n",
       "       [ 0.,  0.,  0.],\n",
       "       [-2.,  0.,  0.]])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Not commutative\n",
    "E32 @ E31 @ E21  - E21 @ E31 @ E32"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c8ce03a5",
   "metadata": {},
   "source": [
    "- $E_{32} E_{31} E_{21} A = U$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "13483e20",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Undo Gauss Eliminate\n",
      "[[ 2.  1.  1.]\n",
      " [ 4. -6.  0.]\n",
      " [-2.  7.  2.]]\n",
      "Lower Triangular matrix\n",
      "[[ 1.  0.  0.]\n",
      " [ 2.  1.  0.]\n",
      " [-1. -1.  1.]]\n"
     ]
    }
   ],
   "source": [
    "iE32 = ematrix(3, 1, 2, -1)\n",
    "iE31 = ematrix(3, 0, 2, -1)\n",
    "iE21 = ematrix(3, 0, 1, 2)\n",
    "\n",
    "iE = iE21 @ iE31 @ iE32\n",
    "print(\"Undo Gauss Eliminate\")\n",
    "print(iE @ U)\n",
    "\n",
    "print(\"Lower Triangular matrix\")\n",
    "print(iE)\n",
    "\n",
    "L =iE"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6c9f9a20",
   "metadata": {},
   "source": [
    "- LU decomposition\n",
    "    * Formula\n",
    "    \n",
    "        $$\n",
    "        LU = \n",
    "        \\begin{bmatrix}\n",
    "        1 & 0 & 0 \\\\\n",
    "        l_{21} & 1 & 0 \\\\\n",
    "        l_{31} & l_{23} & 1\\\\\n",
    "        \\end{bmatrix}\n",
    "        \\begin{bmatrix}\n",
    "        u_{11} & u_{12} & u_{13} \\\\\n",
    "        0 & u_{22} & u_{23} \\\\\n",
    "        0 & 0 & u_{33} \\\\\n",
    "        \\end{bmatrix}\n",
    "        $$\n",
    "        \n",
    "        - $l_{ij} = a'_{ij} / a'_{ii}$\n",
    "        \n",
    "    * Save in a matrix\n",
    "        \n",
    "        $$\n",
    "        M = \\begin{bmatrix}\n",
    "        u_{11} & u_{12} & u_{13} \\\\\n",
    "        l_{21} & u_{22} & u_{23} \\\\\n",
    "        l_{31} & l_{23} & u_{33} \\\\\n",
    "        \\end{bmatrix}\n",
    "        $$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9e15cfd6",
   "metadata": {},
   "source": [
    "## Implementation\n",
    "- Make $A \\rightarrow M$\n",
    "    \n",
    "    $$\n",
    "    U = \\left[\\begin{array}{ccc}\n",
    "    2 & 1 & 1\\\\\n",
    "    0 & -8 & -2\\\\\n",
    "    0 & 0 & 1\n",
    "    \\end{array}\\right],\n",
    "    L = \\left[\\begin{array}{ccc}\n",
    "    1 & 0 & 0\\\\\n",
    "    2 & 1 & 0\\\\\n",
    "    -1 & -1 & 1\n",
    "    \\end{array}\\right]\n",
    "    $$\n",
    "    \n",
    "    $$\n",
    "    M = \\left[\\begin{array}{ccc}\n",
    "    2 & 1 & 1\\\\\n",
    "    2 & -8 & -2\\\\\n",
    "    -1 & -1 & 1\n",
    "    \\end{array}\\right],\n",
    "    $$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "01c8e418",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[2, 1, 1], [4, -6, 0], [-2, 7, 2]])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "9f03bdac",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.0 [ 2 -8 -2]\n"
     ]
    }
   ],
   "source": [
    "# First pivot a[0, 0]\n",
    "# eliminate a[1, 0]\n",
    "i = 0\n",
    "j = 1\n",
    "ratio = A[j, i] / A[i, i]\n",
    "\n",
    "A[j, i] = ratio\n",
    "A[j, i+1:] = A[j, i+1:] - ratio*A[i, i+1:]\n",
    "print(ratio, A[j])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "d966d262",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-1.0 [-1  8  3]\n"
     ]
    }
   ],
   "source": [
    "# eliminate a[2, 0]\n",
    "j = 2\n",
    "ratio = A[j, i] / A[i, i]\n",
    "\n",
    "A[j, i] = ratio\n",
    "A[j, i+1:] = A[j, i+1:] - ratio*A[i, i+1:]\n",
    "print(ratio, A[j])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "64d8ce3c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-1.0 [-1 -1  1]\n"
     ]
    }
   ],
   "source": [
    "# next pivot a_{2,2}\n",
    "# eliminate a_{3, 2}\n",
    "i = 1\n",
    "j = 2\n",
    "\n",
    "ratio = A[j, i] / A[i, i]\n",
    "\n",
    "A[j, i] = ratio\n",
    "A[j, i+1:] = A[j, i+1:] - ratio*A[i, i+1:]\n",
    "print(ratio, A[j])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "bb1ba558",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2  1  1]\n",
      " [ 2 -8 -2]\n",
      " [-1 -1  1]]\n"
     ]
    }
   ],
   "source": [
    "print(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e1a3511",
   "metadata": {},
   "source": [
    "- Substibution\n",
    "    * Forward : $Lc=b$\n",
    "    * Backward: $Ux=c$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "f64d293f",
   "metadata": {},
   "outputs": [],
   "source": [
    "b = np.array([5, -2, 9])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "13ddb204",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "np.int64(5)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Forward (b to c)\n",
    "# First row\n",
    "b[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "8fa59f5c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Second row\n",
    "b[1] = b[1] - A[1, 0]*b[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "319babde",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[  5 -12   2]\n"
     ]
    }
   ],
   "source": [
    "# Third row\n",
    "i = 2\n",
    "\n",
    "for j in range(i):\n",
    "    b[i] -= A[i, j]*b[j]\n",
    "    \n",
    "print(b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "f53ab74d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Back substitution (Same as Gauss Elimination)\n",
    "x = np.empty_like(b)\n",
    "\n",
    "# Third row\n",
    "i = 2\n",
    "x[i] = b[i] / A[i,i]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "16509053",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Second row\n",
    "i = 1\n",
    "x[i] = (b[i] - A[i, i+1]*x[i+1]) / A[i, i]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "1a8a1a06",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1 1 2]\n"
     ]
    }
   ],
   "source": [
    "# First row\n",
    "n = 3\n",
    "i = 0\n",
    "x[i] = b[i]\n",
    "\n",
    "for j in range(i+1, n):\n",
    "    x[i] -= A[i, j]*x[j]\n",
    "    \n",
    "x[i] /= A[i,i]\n",
    "\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "024c240c",
   "metadata": {},
   "source": [
    "### DIY\n",
    "다음 2개의 함수를 만드시오\n",
    "* LU 분해: LUdecomp\n",
    "* Substitution : LUsubs"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9700a890",
   "metadata": {},
   "source": [
    "## 역행렬\n",
    "역행렬은 $Ax_i = e_i$ 를 푸는 과정이다.\n",
    "\n",
    "$$\n",
    "A X = I = A \n",
    "\\begin{bmatrix}\n",
    "x_1 & x_2 & ... & x_n\n",
    "\\end{bmatrix} =\n",
    "\\begin{bmatrix}\n",
    "e_1 & e_2 & ... & e_n\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "### 예제\n",
    "다음 $A$ 행렬의 역행렬을 구하시오.\n",
    "\n",
    "$$\n",
    "A = \\begin{bmatrix}\n",
    "2 & 1 & 1 \\\\\n",
    "4 & 6 & 0 \\\\\n",
    "-2 & 7 & 2\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "        \n",
    "LU 분해 결과는 아래와 같다.\n",
    "\n",
    "$$\n",
    "M = \\left[\\begin{array}{ccc}\n",
    "2 & 1 & 1\\\\\n",
    "2 & -8 & -2\\\\\n",
    "-1 & -1 & 1\n",
    "\\end{array}\\right],\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "1ed5aa61",
   "metadata": {},
   "outputs": [],
   "source": [
    "# LU 분해 결과 (이전 계산)\n",
    "A = np.array([[ 2,  1,  1],  [ 2, -8, -2],  [-1, -1,  1]], dtype=np.float64)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "c7c076ba",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Allocate Identity and inverse Matrix\n",
    "I = np.eye(3)\n",
    "Ainv = np.empty_like(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "72e9f0df",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0.75  0.5  -1.  ]\n"
     ]
    }
   ],
   "source": [
    "# For first column\n",
    "b = I[:, 0].copy()\n",
    "x = Ainv[:, 0]\n",
    "\n",
    "# Forward (b to c)\n",
    "# Second row\n",
    "b[1] = b[1] - A[1, 0]*b[0]\n",
    "\n",
    "# Third row\n",
    "i = 2\n",
    "\n",
    "for j in range(i):\n",
    "    b[i] -= A[i, j]*b[j]\n",
    "    \n",
    "# Back substitution (Same as Gauss Elimination)\n",
    "# Third row\n",
    "i = 2\n",
    "x[i] = b[i] / A[i,i]\n",
    "\n",
    "# Second row\n",
    "i = 1\n",
    "x[i] = (b[i] - A[i, i+1]*x[i+1]) / A[i, i]\n",
    "\n",
    "# First row\n",
    "n = 3\n",
    "i = 0\n",
    "x[i] = b[i]\n",
    "\n",
    "for j in range(i+1, n):\n",
    "    x[i] -= A[i, j]*x[j]\n",
    "    \n",
    "x[i] /= A[i,i]\n",
    "\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "1eb645c5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-0.3125 -0.375   1.    ]\n"
     ]
    }
   ],
   "source": [
    "# For first column\n",
    "b = I[:, 1].copy()\n",
    "x = Ainv[:, 1]\n",
    "\n",
    "# Forward (b to c)\n",
    "# Second row\n",
    "b[1] = b[1] - A[1, 0]*b[0]\n",
    "\n",
    "# Third row\n",
    "i = 2\n",
    "\n",
    "for j in range(i):\n",
    "    b[i] -= A[i, j]*b[j]\n",
    "    \n",
    "# Back substitution (Same as Gauss Elimination)\n",
    "# Third row\n",
    "i = 2\n",
    "x[i] = b[i] / A[i,i]\n",
    "\n",
    "# Second row\n",
    "i = 1\n",
    "x[i] = (b[i] - A[i, i+1]*x[i+1]) / A[i, i]\n",
    "\n",
    "# First row\n",
    "n = 3\n",
    "i = 0\n",
    "x[i] = b[i]\n",
    "\n",
    "for j in range(i+1, n):\n",
    "    x[i] -= A[i, j]*x[j]\n",
    "    \n",
    "x[i] /= A[i,i]\n",
    "\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "940f54b2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-0.375 -0.25   1.   ]\n"
     ]
    }
   ],
   "source": [
    "# For first column\n",
    "b = I[:, 2].copy()\n",
    "x = Ainv[:, 2]\n",
    "\n",
    "# Forward (b to c)\n",
    "# Second row\n",
    "b[1] = b[1] - A[1, 0]*b[0]\n",
    "\n",
    "# Third row\n",
    "i = 2\n",
    "\n",
    "for j in range(i):\n",
    "    b[i] -= A[i, j]*b[j]\n",
    "    \n",
    "# Back substitution (Same as Gauss Elimination)\n",
    "# Third row\n",
    "i = 2\n",
    "x[i] = b[i] / A[i,i]\n",
    "\n",
    "# Second row\n",
    "i = 1\n",
    "x[i] = (b[i] - A[i, i+1]*x[i+1]) / A[i, i]\n",
    "\n",
    "# First row\n",
    "n = 3\n",
    "i = 0\n",
    "x[i] = b[i]\n",
    "\n",
    "for j in range(i+1, n):\n",
    "    x[i] -= A[i, j]*x[j]\n",
    "    \n",
    "x[i] /= A[i,i]\n",
    "\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "8ea60888",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.75   -0.3125 -0.375 ]\n",
      " [ 0.5    -0.375  -0.25  ]\n",
      " [-1.      1.      1.    ]]\n"
     ]
    }
   ],
   "source": [
    "print(Ainv)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "75c70c63",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1., 0., 0.],\n",
       "       [0., 1., 0.],\n",
       "       [0., 0., 1.]])"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Validate\n",
    "A = np.array([[2,1,1],[4,-6,0],[-2,7,2]])\n",
    "A @ Ainv"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "895b7896",
   "metadata": {},
   "source": [
    "### DIY\n",
    "앞서 만든 LU 분해 함수를 이용하여 역행렬 계산 함수를 만드시오"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "84e067bf",
   "metadata": {},
   "source": [
    "## 특수 행렬 분해\n",
    "### Cholesky Decomposition\n",
    "- $LU$ 분해에서 $U$ 행렬은 Diagonal 항과 Off-diagonal 항으로 분리할 수 있다.\n",
    "\n",
    "    $$\n",
    "        A = LU = LDU\n",
    "    $$\n",
    "    \n",
    "- $A$ 행렬이 대칭 (symmetric) 일 경우 다음 관계를 만족한다.\n",
    "\n",
    "    $$\n",
    "        A = A^T \\\\\n",
    "        LDU = (LDU)^T = U^T D L^T\n",
    "    $$\n",
    "    \n",
    "    * 즉 $L = U^T$ 를 만족한다.\n",
    "    * $A=LDL^T$\n",
    "    \n",
    "- Cholesky 분해법\n",
    "    * A 가 Positive Definite 일 때 (symmetric, all positive pivots)\n",
    "    \n",
    "    $$\n",
    "        A = LDL^T = (LD^{1/2}) (D^{1/2} L^T) = (LD^{1/2}) (LD^{1/2})^T\n",
    "    $$\n",
    "    \n",
    "- 구현\n",
    "    \n",
    "    $$\n",
    "    \\begin{align}\n",
    "    l_{ki} &= \\frac{a_{ki} - \\sum_{j=1}^{i-1} l_{ij} l_{kj}}{l_{ii}} \\textrm{ for } i=1,2,...k-1 \\\\\n",
    "    l_{kk} &= \\sqrt{a_{kk} - \\sum_{j=1}^{k-1} l_{kj}^2}\n",
    "    \\end{align}\n",
    "    $$\n",
    "    \n",
    "#### 예제\n",
    "다음 대칭 행렬을 Cholesky 분해하시오.\n",
    "\n",
    "$$\n",
    "A = \\begin{bmatrix}\n",
    "6 & 15 & 55 \\\\\n",
    "15 & 55 & 225 \\\\\n",
    "55 & 225 & 979\n",
    "\\end{bmatrix}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "dab369ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[6, 15, 55], [15, 55, 225], [55, 225, 979]], dtype=np.float64)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "244f0142",
   "metadata": {},
   "outputs": [],
   "source": [
    "L = np.zeros_like(A)\n",
    "\n",
    "# first row\n",
    "k = 0\n",
    "L[k, k] = np.sqrt(A[k, k])\n",
    "\n",
    "# Second row\n",
    "k = 1\n",
    "i = 0\n",
    "L[k, i] = A[k, i] / L[i, i]\n",
    "L[k, k] = np.sqrt(A[k, k] - L[k, i]**2)\n",
    "\n",
    "# Third row\n",
    "k = 2\n",
    "i = 0 \n",
    "L[k, i] = A[k, i] / L[i, i]\n",
    "\n",
    "i = 1\n",
    "L[k, i] = (A[k, i] - sum(L[i, j]*L[k, j] for j in range(i)))/L[i,i]\n",
    "\n",
    "L[k ,k] = np.sqrt(A[k, k] - sum(L[k, j]**2 for j in range(k)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "7eebca22",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 2.44948974,  0.        ,  0.        ],\n",
       "       [ 6.12372436,  4.18330013,  0.        ],\n",
       "       [22.45365598, 20.91650066,  6.11010093]])"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "L"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "032d669a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[  6.,  15.,  55.],\n",
       "       [ 15.,  55., 225.],\n",
       "       [ 55., 225., 979.]])"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Validate\n",
    "L @ L.T"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6efc5570",
   "metadata": {},
   "source": [
    "- LU decomposition과 비교"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "id": "72cf5f99",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 6.         15.         55.        ]\n",
      " [ 2.5        17.5        87.5       ]\n",
      " [ 9.16666667  5.         37.33333333]]\n"
     ]
    }
   ],
   "source": [
    "# LU decomposition\n",
    "M = A.copy()\n",
    "i = 0\n",
    "j = 1\n",
    "ratio = M[j, i] / M[i, i]\n",
    "\n",
    "M[j, i] = ratio\n",
    "M[j, i+1:] = A[j, i+1:] - ratio*M[i, i+1:]\n",
    "\n",
    "j = 2\n",
    "ratio = M[j, i] / M[i, i]\n",
    "\n",
    "M[j, i] = ratio\n",
    "M[j, i+1:] = M[j, i+1:] - ratio*M[i, i+1:]\n",
    "\n",
    "i = 1\n",
    "j = 2\n",
    "ratio = M[j, i] / M[i, i]\n",
    "\n",
    "M[j, i] = ratio\n",
    "M[j, i+1:] = M[j, i+1:] - ratio*M[i, i+1:]\n",
    "\n",
    "print(M)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "112200aa",
   "metadata": {},
   "source": [
    "- $A = LU=LDU'$\n",
    "\n",
    "$$\n",
    "L = \\begin{bmatrix}\n",
    "1 & 0 & 0 \\\\\n",
    "2.5 & 1 & 0 \\\\\n",
    "9.1667 & 5 & 1\n",
    "\\end{bmatrix}\n",
    "U = \\begin{bmatrix}\n",
    "6 & 1.5 & 55 \\\\\n",
    "0 & 17.5 & 87.5 \\\\\n",
    "0 & 0 & 37.3333\n",
    "\\end{bmatrix}\n",
    "= \\begin{bmatrix}\n",
    "6 & 0 & 0 \\\\\n",
    "0 & 17.5 & 0 \\\\\n",
    "0 & 0& 37.3333\n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "1 & 2.5 & 9.1667 \\\\\n",
    "0 & 1 &  5 \\\\\n",
    "0 & 0 & 1\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "- $L\\sqrt{D}$\n",
    "\n",
    "$$\n",
    "L \\sqrt{D} = \\begin{bmatrix}\n",
    "\\sqrt{6} & 0 & 0 \\\\\n",
    "2.5 \\sqrt{6} & \\sqrt{17.5} & 0 \\\\\n",
    "9.1667 \\sqrt{6} & 5 \\sqrt{17.5} & \\sqrt{37.3333}\n",
    "\\end{bmatrix}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c509bc0",
   "metadata": {},
   "source": [
    "### 삼중대각 행렬\n",
    "- 삼중 대각 행렬 (Tri-diagonal matrix)는 Diagonal과 그 바로 옆에 1개씩만 값이 있는, 즉 3의 대역폭을 갖는 행렬이다.\n",
    "\n",
    "- 실제 편미분 방정식 해석에서 자주 보이는 행렬이다.\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix}\n",
    "b_1 & c_1 & & & & & \\\\\n",
    "a_2 & b_2 & c_2 & & & & \\\\\n",
    "& a_3 & b_3 & c_3 & & & \\\\\n",
    "& & & & & & & \\\\\n",
    "& & & & a_{n-1} & b_{n-1} & c_{n-1} \\\\\n",
    "& & & & & a_{n} & b_{n}\n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "x_1 \\\\ x_2 \\\\ x_3 \\\\ \\\\ x_{n-1} \\\\ x_n\n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "d_1 \\\\ d_2 \\\\ d_3 \\\\ \\\\ d_{n-1} \\\\ d_n\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "- (LU Decomposition + Forward substitution), Backward substitution 으로 효과적으로 계산할 수 있다.\n",
    "\n",
    "- Algorithm\n",
    "    * For $i = 2, 3,...,n$\n",
    "        \n",
    "        $$\n",
    "        \\begin{align}\n",
    "            w &= \\frac{a_i}{b_{i-1}} \\\\\n",
    "            b_i &= b_i - w c_{i-1} \\\\\n",
    "            d_i &= d_i - w d_{i-1}\n",
    "         \\end{align}\n",
    "        $$\n",
    "        \n",
    "    * Back substitution\n",
    "    \n",
    "        $$\n",
    "        \\begin{align}\n",
    "        x_n &= \\frac{d_n}{b_n} \\\\\n",
    "        x_i & = \\frac{d_i - c_i x_{i+1}}{b_i} \\textrm{ for } i=n-1, n-2,...,1        \n",
    "        \\end{align}\n",
    "        $$\n",
    "\n",
    "#### 예제\n",
    "아래 삼중대각 행렬의 해를 구하시오.\n",
    "\n",
    "$$\n",
    "\\begin{bmatrix}\n",
    "2.04 & -1 & & \\\\\n",
    "-1 & 2.04 & -1 & \\\\\n",
    "& -1 & 2.04 & -1 \\\\\n",
    "& & -1 & 2.04\n",
    "\\end{bmatrix}\n",
    "T \n",
    "= \n",
    "\\begin{bmatrix}\n",
    "40.8 \\\\ 0.8 \\\\ 0.8 \\\\ 200.8\n",
    "\\end{bmatrix}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "id": "f7df4b3e",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[2.04, -1, 0, 0], [-1, 2.04, -1, 0], [0, -1, 2.04, -1], [0, 0, -1, 2.04]])\n",
    "d = np.array([40.8, 0.8, 0.8, 200.8])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "id": "fa9a4e0c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 2.04, -1.  ,  0.  ,  0.  ],\n",
       "       [-1.  ,  2.04, -1.  ,  0.  ],\n",
       "       [ 0.  , -1.  ,  2.04, -1.  ],\n",
       "       [ 0.  ,  0.  , -1.  ,  2.04]])"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "1f0ebcf3",
   "metadata": {},
   "outputs": [],
   "source": [
    "a = np.array([0., -1., -1. , -1.])\n",
    "b = np.array([2.04, 2.04, 2.04, 2.04])\n",
    "c = np.array([-1., -1., -1., 0.])\n",
    "\n",
    "n = len(a)\n",
    "x = np.empty_like(a)\n",
    "\n",
    "# Forward sweep\n",
    "for i in range(1, n):\n",
    "    w = a[i] / b[i-1]\n",
    "    b[i] -= w*c[i-1]\n",
    "    d[i] -= w*d[i-1]\n",
    "\n",
    "# Backward sweep    \n",
    "x[-1] = d[-1] / b[-1]\n",
    "for i in range(n-2, -1, -1):\n",
    "    x[i] = (d[i] - c[i]*x[i+1]) / b[i]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "id": "a5150aa6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 65.96983437,  93.77846211, 124.53822833, 159.47952369])"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Solution\n",
    "x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "id": "5db24d68",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 40.8,   0.8,   0.8, 200.8])"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Validation\n",
    "A @ x"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a39720d1",
   "metadata": {},
   "source": [
    "### DIY\n",
    "- Cholesky 분해 계산 함수를 만드시오.\n",
    "- 삼중 대각 행렬 계산 함수를 만드시오. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bd25f5b8",
   "metadata": {},
   "source": [
    "## Scipy 활용\n",
    "\n",
    "`scipy.linalg` 모듈은 다양한 행렬 계산 알고리즘을 제공함\n",
    "* `numpy.linalg` 는 제한된 기능을 제공\n",
    "\n",
    "**참고**\n",
    "* https://scipy-lectures.org/intro/scipy.html#linear-algebra-operations-scipy-linalg\n",
    "* https://docs.scipy.org/doc/scipy/tutorial/linalg.html\n",
    "\n",
    "### Solve linear system\n",
    "\n",
    "- `linalg.solve` : $Ax=b$ 해석\n",
    "- `linalg.inv` : $A^{-1}$ 계산\n",
    "- `linalg.det` : $det(A)$ 계산\n",
    "- `linalg.norm` : 벡터, 행렬의 Norm 계산\n",
    "\n",
    "### Decomposition\n",
    "- LU (Lower Upper) Decomposition\n",
    "- Singual Value Decomposition (SVD)\n",
    "- QR Decomposition\n",
    "\n",
    "### Eigenvalue \n",
    "- `linalg.eig` : Eigenvalue, Eigenvector 계산\n",
    "    * $Ax = \\lambda x$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "id": "fbecfb26",
   "metadata": {},
   "outputs": [],
   "source": [
    "np.linalg.solve?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "id": "7564722b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[ 4. , -6. ,  0. ],\n",
       "        [ 0.5,  4. ,  1. ],\n",
       "        [-0.5,  1. ,  1. ]]),\n",
       " array([1, 1, 2], dtype=int32))"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# LU decompse\n",
    "A = np.array([[2,1,1],[4,-6,0],[-2,7,2]])\n",
    "\n",
    "from scipy import linalg\n",
    "linalg.lu_factor(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "id": "a8c46fd4",
   "metadata": {},
   "outputs": [],
   "source": [
    "linalg.lu_factor?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3c1f5088",
   "metadata": {},
   "source": [
    "### 예제\n",
    "양 끝단이 고정된 보의 방정식은 다음과 같다.\n",
    "\n",
    "$$\n",
    "u''(x) = 1, u(0) = 0, u(1) = 0.\n",
    "$$\n",
    "\n",
    "이 미분방정식의 이론해는 다음과 같다.\n",
    "\n",
    "$$\n",
    "u(x) = -\\frac{1}{2}x + \\frac{1}{2}x^2.\n",
    "$$\n",
    "\n",
    "수치적으로 이를 해석하는 경우 $[0,1]$ 구간을 *n* 등분하고, 각 지점의 값을 $u_i = u(x_i)$ 라 하면 위 미분 방정식으 다음 형태로 바뀐다.\n",
    "\n",
    "$$\\frac{1}{\\Delta x^2}\n",
    "\\begin{bmatrix}\n",
    "-2 & 1 & 0 &... &0 \\\\\n",
    "1 & -2& 1 & 0 ... & 0 \\\\\n",
    "0 & 1 & -2 & 1 ... & 0 \\\\\n",
    "... & ... & ... & ... & ... \\\\\n",
    "0 & 0 & 1 & -2 & 1 \\\\\n",
    "0 & 0 & 0 & 1 & -2\n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "u_1 \\\\\n",
    "u_2 \\\\\n",
    "u_3 \\\\\n",
    "... \\\\\n",
    "u_{n-2} \\\\\n",
    "u_{n-1}\n",
    "\\end{bmatrix}\n",
    "= \\begin{bmatrix}\n",
    "f(x_1) \\\\\n",
    "f(x_2) \\\\\n",
    "f(x_3) \\\\\n",
    "... \\\\\n",
    "f(x_{n-2})\\\\\n",
    "f(x_{n-1})\\\\\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "이 문제를 수치적으로 해석해보자."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "id": "c29bcc98",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Number of division\n",
    "n = 51\n",
    "x = np.linspace(0, 1, n+1)\n",
    "h = np.diff(x)[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "id": "ba41eba7",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Matrix system\n",
    "A = np.diag([-2]*(n-1)) + np.diag([1]*(n-2), -1) + np.diag([1]*(n-2), 1)\n",
    "b = np.array([1]*(n-1))*h**2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "8e6592d7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[-2,  0,  0, ...,  0,  0,  0],\n",
       "        [ 0, -2,  0, ...,  0,  0,  0],\n",
       "        [ 0,  0, -2, ...,  0,  0,  0],\n",
       "        ...,\n",
       "        [ 0,  0,  0, ..., -2,  0,  0],\n",
       "        [ 0,  0,  0, ...,  0, -2,  0],\n",
       "        [ 0,  0,  0, ...,  0,  0, -2]], shape=(50, 50)),\n",
       " array([[0, 0, 0, ..., 0, 0, 0],\n",
       "        [1, 0, 0, ..., 0, 0, 0],\n",
       "        [0, 1, 0, ..., 0, 0, 0],\n",
       "        ...,\n",
       "        [0, 0, 0, ..., 0, 0, 0],\n",
       "        [0, 0, 0, ..., 1, 0, 0],\n",
       "        [0, 0, 0, ..., 0, 1, 0]], shape=(50, 50)),\n",
       " array([[0, 1, 0, ..., 0, 0, 0],\n",
       "        [0, 0, 1, ..., 0, 0, 0],\n",
       "        [0, 0, 0, ..., 0, 0, 0],\n",
       "        ...,\n",
       "        [0, 0, 0, ..., 0, 1, 0],\n",
       "        [0, 0, 0, ..., 0, 0, 1],\n",
       "        [0, 0, 0, ..., 0, 0, 0]], shape=(50, 50)))"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.diag([-2]*(n-1)), np.diag([1]*(n-2), -1), np.diag([1]*(n-2), 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "id": "ed8279c2",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solve\n",
    "#u = np.linalg.solve(A, b)\n",
    "u = linalg.solve(A, b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "id": "41d10132",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0xffd91179dbd0>"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1QAAAJqCAYAAAAsZgE0AAAAQHRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjErZGZzZzEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvzRIYmAAAAAlwSFlzAAAXEgAAFxIBZ5/SUgAAo05JREFUeJzs3Xd0HNXdPvBne1HvXbJlWW6SiyxbcpdxxTamhhDiGJKQ5CUvAQeSkJCAf5AChGB4SeJUErCpBmzce5O7bLnJTa6ybPXeV1t/fyyalVCx+t1dPZ9zOGfu7MzoWXss9rv3zr0ym81mAxEREREREXWZXHQAIiIiIiIiV8WCioiIiIiIqJtYUBEREREREXUTCyoiIiIiIqJuYkFFRERERETUTSyoiIiIiIiIuokFFRERERERUTexoCIiIiIiIuomFlRERERERETdxIKKiIiIiIiom1hQERERERERdRMLKiIiIiIiom5iQUVERERERNRNLKiIiIiIiIi6SSk6wEBVWFgIm80m7OcHBgYCAEpLS4VlIOfF+4Paw3uD2sN7g9rDe4M6Ivr+kMlkCA0N7dE1WFAJYrPZhBZUzXMQtYf3B7WH9wa1h/cGtYf3BnXEle8PDvkjIiIiIiLqJhZURERERERE3cSCioiIiIiIqJtYUBEREREREXUTCyoiIiIiIqJuYkFFRERERETUTSyoiIiIiIiIuokFFRERERERUTexoCIiIiIiIuomFlRERERERETdxIKKiIiIiIiom1hQERERERERdRMLKiIiIiIiom5iQUVERERERNRNLKiIiIiIiIi6SSk6wNcZjUZ8+eWXOHToEEpLS+Hp6YkxY8bgm9/8JgICArp0rbq6Onz22WfIyMhAZWUlfH19MWHCBDz88MPw8PBo8xyr1YqtW7diz549KCwshFarxahRo/Dwww8jMjKyN94iERERERG5CafqoTIajfjtb3+Lzz//HAaDAcnJyQgICMC+ffvw/PPPo7CwsNPXqqmpwQsvvIAtW7ZAoVBgwoQJ0Ol02Lp1K371q1+hpqam1Tk2mw1vv/023n//fZSXlyMpKQlRUVE4duwYfvnLX+LKlSu9+XaJiIiIiMjFOVUP1bp165CdnY34+Hj85je/gVarBQBs2rQJq1atwt/+9je8/PLLnbrW+++/j4KCAkycOBE//elPoVAoAAD/+c9/sG3bNrz//vt46qmnWpyzd+9eHD16FGFhYXj55Zfh6+sLADh69ChWrFiBd955B2+//bZ0LSIiIiIiGticpofKbDZj27ZtAIDvf//7UjEFAIsWLUJMTAwuXryI69ev3/FalZWVOHDgABQKBZ544okWBdB3vvMdeHt74+DBg6isrGxx3qZNmwAA3/72t6ViCgBSU1ORnJyMoqIiHD9+vAfvkoiIiIiI3InTFFSXLl1CXV0dQkJCMHjw4Favp6SkAABOnDhxx2udOnUKNpsNI0eObFEYAYBKpcL48eNhtVpx+vRpaX9xcTFu374NtVqNpKSkVtdMTU0FAGRmZnbhXRERERERkTtzmoLq5s2bANBmMQUAsbGxLY7rybWa9ufk5Ej7mrajoqKgVLYeCdl0Tmd+vrOz2Wy4VlqHfx7OQWZereg4RERERDQApedU47/HcnGrokF0lB5xmmeoSktLAaDdmfz8/f1bHNeZazWd83VNP6P5te7089s6pyPPPvtsq31qtRqvvfYaACAwMLBT1+kLfz+Ug1XHswEA04cEYP7YtgtPGriavlQICgoSnIScDe8Nag/vDWoP7w1qz/Y9eThXUIN/HbmJp6fH4pGkCNGRusVpeqgMBgMAQKPRtPl60zNVTcf15FpN+xsbG1udo1arOzynMz/f2Y2N8JG2j+aUo67RLDANEREREQ00hdUGnCtwzLo9LtKng6Odm9P0UNlsth693taxMpmsyzm6c05bVqxY0eHrpaWlXXpPvSlGZ4OPVokqgxlGiw1bzuQgbbDr3sTU+5q+RSwpKRGchJwN7w1qD+8Nag/vDWrLxotl0nakrxZ+qEdJSf8P/ZPJZAgLC+vRNZymh0qn0wFo2WvUXNP+5rP/3ela7fUmNV2reQ9W03V74+c7O6VchhlxjiGHB29WC0xDRERERAPNwZuO3qlZ8UG91qkhgtMUVE3PFJWVlbX5enl5eYvjOnOtpnO+rulnNL/WnX5+W+e4stnxjvdxqqAOtY0WgWmIiIiIaKAorDHiSpmj42N2vGs/X+c0BVVMTAwA4MaNG22+3rT+VNNxPblW0/7m1xo0aBAA4NatWzCbWz9T1HROdHT0HX++Kxgb6Qs/vQoAYLYCR2/X3OEMIiIiIqKeO5Tr+Nw5yF+H2AC9wDQ95zQF1fDhw6HX61FUVNRmIXTs2DEAaHONqK8bO3YsZDIZLl68iKqqqhavmUwmZGZmQiaTYdy4cdL+4OBgREREwGg04uTJk62uefToUQDA+PHju/S+nJVSLkNas2F/h26yoCIiIiKivtf8cZO7hrr2cD/AiQoqpVKJ+fPnAwD+85//tHj+adOmTbh58yaGDx+OuLg4af+2bduwbNkyfPTRRy2u5efnhylTpsBsNuPf//43LBbHcLYPPvgA1dXVmDp1aqtFfxctWgQA+PDDD1sUYseOHcOJEycQHByMCRMm9Np7Fq35sL8zhXWo5rA/IiIiIupD+dVGXK9wzFkwe5hrD/cDnGiWPwB44IEHkJWVhezsbDzzzDMYPnw4SktLceXKFXh5eeHHP/5xi+Orq6uRn5+PioqKVtd6/PHHceXKFRw7dgzLli3DkCFDcOvWLdy6dQshISF47LHHWp0zc+ZMnDp1ChkZGVi2bBkSExNRU1ODCxcuQKVS4Sc/+Umbi/66qtHhPvDTKVHRYIbFBhy9VYO5cb6iYxERERGRm2reOxXjq8Egf9ce7gc4UQ8VYF8Davny5XjwwQehVqtx/PhxFBcXY8aMGXj99dcRGhra6Wt5e3vj1Vdfxfz582E2m5GRkYH6+nrMnz8fr776Kry9vVudI5fL8eyzz2Lp0qXw9/dHZmYmcnNzMWHCBLz22msYNmxYb75d4RRyGaZEe0ntA5ztj4iIiIj60MFmz09NjfHq4EjXIbOJWgxpgCsoKBC2DhXgWBMi/fxN/HJnLgBALgP+e38cfHXu0wtH3cM1Q6g9vDeoPbw3qD28N6jJrapGPLXJMVfC3+6JxZghEQDE3R9utQ4ViTEsSIcAvb2AstqAI7c4OQURERER9b7mw/1i/TQI91YLTNN7WFANcHKZDFObDfvjIr9ERERE1NtsNluLxXynxrR+/MZVsaCiFjf0+eIGlNWbBKYhIiIiIndzs7IRt6uNUttdnp8CWFARgKEBWgR72Bf5tQE4nMthf0RERETUe5r3Tg0N0CLE0z2G+wEsqAj2h/Gaf0twiAUVEREREfUSm82Gg7mOx0rcqXcKYEFFX2k+7O9iSQNK6jjsj4iIiIh67kZFIwpqHJ8tp0S7z/NTAAsq+kqsnwZhXiqpzWF/RERERNQbmq91OjxQhyAPVQdHux4WVATAPuyv+bcFXOSXiIiIiHqq9ex+7jXcD2BBRc1Ma3aDXykzoKjW2MHRREREREQdu1JmQPFXj5LIAEyOZkFFbizGV4PIZgusNf82gYiIiIioq5pPdjYyWIcAvXsN9wNYUFEzX5/tj4v8EhEREVF3WW22Fo+RuNtkFE1YUFELzWf7u17RiNvVjQLTEBEREZGruljSgLJ6MwBALgOmuOHzUwALKvqaKB8NBvtppPbBHA77IyIiIqKuO5Dj6J0aHaKHr1YpME3fYUFFrUxr1kuVfrMaNptNYBoiIiIicjUWq63F81PTBrnncD+ABRW1oflzVHnVRtyo4LA/IiIiIuq8M4V1qG60AACUchlSo9xzuB/AgoraEOKpxrBAndTmmlRERERE1BUHms0WPT7cA55qhcA0fYsFFbWp+ZpUB3I47I+IiIiIOsdoseLorWbD/WLcd7gfwIKK2jE1xhtymX27pN6MS6UNYgMRERERkUs4mV+HepMVAKBRyDAh0lNwor7Fgora5KdTIiFYL7UPcJFfIiIiIuqE9Gaz+6VEekGrdO+Sw73fHfVI89lYDt2shsXKYX9ERERE1L4GkxXH82ql9tRB7jsZRRMWVNSuSVFeaPpCodJgwbnierGBiIiIiMipZdyugdFi/xLeQy1HUpiH4ER9jwUVtctLo8C4Zv8ImnffEhERERF9XfPZoSdFeUGlcP9yw/3fIfVI81lZjtyqgcnCYX9ERERE1FpNowWnCuqktrvP7teEBRV1aGKkF9QK+3R/dUYrThXU3uEMIiIiIhqIjtyqgdk+uR98tQokhug7PsFNsKCiDulUckyIcEx1ydn+iIiIiKgtzYf7TYn2gqJpDR43x4KK7qj5bH8Zt2vQ2PTVAxERERERgIoGM84VOSYwa/750d2xoKI7Gh/uAb3KfqsYzDZk3GYvFRERERE5HLxZjaYVdoL0SgwL1IkN1I9YUNEdqRVypEQ0W+Q3/ZTANERERETkbA4cPS9tT432hFw2MIb7ASyoqJOmy0qk7UxZIGqrOIU6EREREQGFtwqQrQqS2tM0A+tzIgsq6pTE5AR4m+zTYJrlShw5nCU4ERERERE5g4PHL0nbEY3lGDwyTmCa/seCijpFpVZhsqpSaqcXGMWFISIiIiKnYLPZkF6pktrTvBshlw+sEmNgvVvqkRmjIqTtc+oQlBUUC0xDRERERKLlXLqGm5pAqT19/MDqnQJYUFEXDBsdjyBjFQDAKpPj4LELghMRERERkUjpp29K23GNJYgYHCUwjRgsqKjTFHI5pns41hdILxs4s7cQERERUUsWixkHDJ5Se3pQBwe7MRZU1CXTk2Kl7avaEORdyREXhoiIiIiEyT55ASVqHwCA3GbF1JSRghOJwYKKumRQXAxijGVSO/3kVYFpiIiIiEiU/ReLpO0EUzECggMEphGHBRV12XR/q7S9v04Pq8UiMA0RERER9TdTYyMOWR0F1IxIncA0YrGgoi6bNnGEtF2g8ce1s5c6OJqIiIiI3M3po2dRo9IDAFRWE1InJQpOJA4LKuqykLBAjDQ6pkzffz5PYBoiIiIi6m/7r1dI28m2Unh66gWmEYsFFXXLtDDHAm4HTX4wG7nQLxEREdFA0FBVjQxFqNSeEesnMI14LKioW6akJkBhtT87VaH2wrmMs4ITEREREVF/OHb4DBoVagCA3mxA0sRRghOJxYKKusXH1wvjrCVSe//Vsg6OJiIiIiJ3kZ7fKG1PVlZAo1Z1cLT7Y0FF3TZjkLe0fVQWgsbaOoFpiIiIiKivVRUU4bQ6TGpPHxnWwdEDAwsq6rYJKaOgtdifnapXanHi8GmxgYiIiIioTx06eh4WuQIA4GeqxajR8YITiceCirpNp9UgRVEutdNvsYeKiIiIyJ2lN3vKY5q+DkoFywn+CVCPTB8WLG1nqsNQU1LSwdFERERE5KqKrl7HRV241J4+LlZgGufBgop6ZMy44fA21wMATHIVjhw+JzgREREREfWF9BNXpe0IUyXi4iIFpnEeLKioR1QKOaZoa6T2gWKLwDRERERE1BdsVisO1Gql9jRfM2QymcBEzoMFFfXYjNHR0naWLhylOTcFpiEiIiKi3pZz5gJu6hyPesyYOExgGufCgop6bPjwQQgxVQMAbDI5DmRcFpyIiIiIiHrT/nO3pe2hplKEhwcJTONcWFBRj8lkMkz3MUrtfdUa2Gw2gYmIiIiIqLeYG43Yb/aX2jNCB/ZCvl/Hgop6RVryUGk7RxeMnKxLAtMQERERUW85n3EG5WpvAIDcZsG0SSMFJ3IuLKioV0RGhSDOWCq192flCkxDRERERL1l71XHuqPjLKXw9fESmMb5sKCiXjMjVClt7zf6wWIyCUxDRERERD1lqK7BUXmI1E4b7C0wjXNiQUW9ZlrqSMht9mnTy9XeOJ+RJTgREREREfVExpEzaFDap0vXWRoxceIowYmcDwsq6jV+ft4YaymR2nuvlHZwNBERERE5u323DdL2JEU5tFq1wDTOiQUV9aq0GMeY2iOyYBjqagWmISIiIqLuqiwswSlNmNSeMTxUYBrnxYKKelVKSgK0lkYAQINSi+OHzgpORERERETdcfDoeVhlCgCAv6kGCWO5mG9bWFBRr9LqNJgkL5Pa+243CExDRERERN21r0wmbU/X10GpYOnQFv6pUK9LG+aYCeaUKhSVxWUdHE1EREREzub2lRxc0Tab3W/cYIFpnBsLKup1CeOGw99UAwCwyBU4ePSc4ERERERE1BX7Tl6TtmMayzB4aIzANM6NBRX1OqVSgWm6Gqm9r6SDg4mIiIjIqVgtFqTXe0jtGf4WgWmcHwsq6hNpYx3dwle0Ici7niswDRERERF1VvaZbBSpfQEAMpsV01NGiA3k5FhQUZ8YHD8I0Y3NJqc4cVVgGiIiIiLqrL0X8qXtBGMxgsKCBKZxfkrRAZrLzs7G2rVrcfnyZZjNZkRGRmLevHlIS0vr1vUyMzOxYcMG5OTkAAAGDRqExYsXY/z48a2OLSkpwYkTJ3D69Gnk5eWhvLwcOp0OsbGxmDdvHpKTk3vwzgYemUyGND8TVtXb2/vr9PiW1Qq5nDU8ERERkbMyNhpxyBIgVQlp4RqxgVyA0xRUGRkZWLFiBWw2G0aMGAEvLy+cO3cOK1euxM2bN/HYY4916XpbtmzBe++9B4VCgcTERCiVSpw9exavv/46Hn/8cSxYsKDF8e+88w6ys7OhVqsxdOhQxMXFoaioCGfOnMGZM2ewcOHCLmcY6KZNHI7Ve8thk8lRpPZF9tnLGDF2uOhYRERERNSOk8eyUKu0Pz+ltpgwaXKi4ETOzykKqtraWqxcuRJWqxXPPfccUlJSAACVlZV46aWXsHnzZowfPx4JCQmdul5+fj5Wr14NlUqF5cuXIz4+Xtr/4osvYvXq1Rg3bhzCwhwrPwcGBmL69OmYNm0atFqttP/kyZN44403sHnzZowdOxZjxozpxXfu3oIjQjHKeAnnNPZVtfedu82CioiIiMiJ7b1eCajsBdVEazE8vFlQ3YlTjL/as2cP6uvrkZycLBVTAODr64slS5YAADZt2tTp623ZsgUWiwVz5syRiikACA8Px/333w+LxYKtW7e2OOeZZ57BnDlzWhRTAJCUlISZM2cCAA4dOtTl9zbQpYWrpe2DFn8YDY0C0xARERFRe2rKK3FC0WztqTg/gWlch1MUVJmZmQCA1NTUVq8lJSVBpVIhKysLRqOxU9c7efJku9ebNGlSi5/ZGTEx9nn3KyoqOn0O2U2eMhpqiwkAUKvUI/PIGcGJiIiIiKgtBw9nwSy3D2DzNtVhbMpowYlcg1MUVLm59im1Y2NjW72mVCoRHR0Nk8mE/Pz8Vq9/XV1dHUpLSwHYJ6H4uoCAAHh5eaGkpAT19fWdyldUVAQA8PHx6dTx5ODh5YkUOBai2pNT08HRRERERCTK3iKrtD1dUwWVyimeDnJ6wv+U6uvrUVdXBwDw9/dv8xh/f39cu3YNpaWlbRZJzTUVUx4eHq2G7zUJCAhATU0NSktLER0d3eH16urqkJ6eDgCYMGFCh8c29+yzz7bap1ar8dprrwGwP7MlklJp/6sPCur7aTAXjYvBgbP2XqpMZSjkZisCwkLucBaJ1J/3B7kW3hvUHt4b1B7eG67hxoXLyNY6Pp/dO310v/yducP9IbyHymAwSNsaTdvTMjbtb37sna7X3rW6er1//etfqK6uxtChQzFx4sQ7Hk+tpU6fAD9TLQDAIldg285jghMRERERUXOb089K2zHGcowcw4nEOqtXeqjefPNN3Lp1q0vnPPXUU4iLi+uNH9+CzWYDYF8Hqae+/PJLHD58GJ6ennj66ae7dM0VK1Z0+HppaamUVYSmbwFKSkrucGTvmK6rwXqzJwBgW74Rc4qLe+XviPpGf98f5Dp4b1B7eG9Qe3hvOD+LxYyd1Vrgq7nEZvhbpFFffU30/SGTyVrM/N0dvVJQlZSUdOr5puYaG+2zvTUfltfY2Ai9Xt+pY9uj0+kAdNz71Jnr7du3Dx9//DE0Gg1++ctfIiSEQ9R6Yub4IVh/zP7M2jVtCG5duoboEb1fUBMRERFR11w4fg4lavtcAXKbFWmTRwlO5Fp6paBqei6oO/R6PfR6Perr61FeXt5mQVVeXg6gc88dNR1TV1cHg8HQZtFUVlbW4fWOHz+Ov//971AoFPjZz37WYup16p7BcdEYfOAgbqjtf+Z7T93AYyyoiIiIiITbm10MKO3zCow2FyMgaKTgRK5F+DNUgGNa8uvXr7d6zWw2Izc3FyqVCuHh4Xe8loeHh1Qo5eTktHq9rKwMNTU1CAwMbLN4O3/+PN5++20AwNNPP82FfHvRzCDH7bav0QfmTk6DT0RERER9w1Bdg8MIltp3xXgKTOOanKKgSkpKAgAcPXq01WsnT56EyWRCQkIC1Gp1q9e7er0jR460OKa569ev449//CPMZjP+53/+p811rKj7pk8ZBbnNAgAoV3vjXMbZO5xBRERERH3p6KHTaFDaR3TpLI1ImZQoOJHrcYqCatasWdDpdDhx4gSOHXPMAFdVVYUPPvgAALBo0aJW5y1btgzLli2ThgQ2WbBgAeRyOXbu3InLly9L+wsKCrBu3TrI5XIsWLCgxTn5+fn4wx/+gIaGBjz++ONIS0vrxXdIAODn54MkS7M1qa6Ud3A0EREREfW1vfmOEUOTFRXQatufKZvaJnwdKgDw9PTEk08+ibfeegsrVqzAyJEj4eXlhaysLNTV1eHuu+9GYmLrarlpIgyz2dxif3h4OJYsWYJVq1Zh+fLlGD16NBQKBc6ePQuj0YilS5e2Gj749ttvo7q6Gt7e3rh+/Tr++te/tvp5ERERuO+++3rvjQ9AM2N9ccK+jjOOKkJRX1EJvZ+v0ExEREREA1Fp7m2c1ThmuLsrMUJgGtflFAUVAKSmpuLll1/G2rVrceXKFZjNZkRERGDevHmYOXNml6+3aNEihIaGYuPGjbh48SIAIDY2FosXL0ZycnKr45sWF66ursb+/fvbvObIkSNZUPXQhJRR8Lh+DnVKHRoVahw+dAazF80QHYuIiIhowNl/9CKssigAQLCpGiMSWn9GpjtzmoIKAIYPH44XXnih08evWbOmw9eTk5PbLJ7a0laPFPU+jVqFKeoq7LDap7ffW2jGbMGZiIiIiAYaq8WCvTVa4KsJsdO8G6GQO8XTQC6Hf2rU7+4aEy1tn9NFoOhajrgwRERERAPQ9TMXcEsbJLVnpgwTmMa1saCifjd8xGCEmaqk9v7jlzs4moiIiIh6296sPGl7mKkU4RHBHRxNHWFBRf1OJpMhzc8xkcjeOg9YvzaxCBERERH1DVN9PdJtgVL7rgjO7NcTLKhIiLTUEdJ2vjYAVzKzBKYhIiIiGjhOHj6FapV9AV+V1YwpkxMEJ3JtLKhIiNAQf4wyN1uT6mKRwDREREREA8eenFppe6KsDF4eOoFpXB8LKhJmZrSHtH1AFoLGmhqBaYiIiIjcX1VBIU5oHOtNzRweIjCNe2BBRcJMSR0FjcW+OnedUoeMQ6cEJyIiIiJyb+mHz8Mst6+c5Guqw7hx8YITuT4WVCSMXqfBJGWF1N592ygwDREREZF7s9ls2FPhWIY2zbMOSgXLgZ7inyAJNWtUuLR9RhuO0pxbAtMQERERua+crIu4rnMM8Zs1YajANO6DBRUJNSpxCIJN1QAAq0yOvRnZghMRERERuafdZxxfXA81lSI6JkxgGvfBgoqEUsjluMunUWrvrdHBarEITERERETkfkwGA/Zbmq09FaYSmMa9sKAi4WY2W5MqTxuA7JPnBaYhIiIicj8nDp9Gtco+w7LKasI0rj3Va1hQkXChYYFIMDnWodpzoVBgGiIiIiL303ztqRRbKby8PDo4mrqCBRU5hbui9NL2QQTDUFvbwdFERERE1FkVRSXIVDmel5o1PEhgGvfDgoqcwuTJidBa7M9S1Su1OHrorOBERERERO5h/+HzsMgVAAB/Uw1GJ424wxnUFSyoyCnodFpMkZdJ7d15jR0cTURERESdYbVaW6w9NVNfB6VSITCR+2FBRU6j+ZpUWeoQFOfmC0xDRERE5Pqun7+Km5pms/tNiBOYxj2xoCKnMWJMPMKMlQAAG9ekIiIiIuqxPWdzpe1hjUWIHBwpMI17YkFFTkMul2Omt0Fq76nWwWq1CkxERERE5LqMBiP2m/2l9l1hyg6Opu5iQUVOZWbqCMhs9iKqUOOLC6cvCU5ERERE5JqOHz2DWqV9JmW1xYipU8YITuSeWFCRUwmOCMFoY7M1qc7zOSoiIiKi7th9w7EMTaqtBJ7engLTuC8WVOR07orSSduHbEFo4JpURERERF1SVliMU6pgqT1reGAHR1NPsKAip5M6dQz0ZvuzVAaFBocOnBGciIiIiMi17Dt8HlaZfXr0QGM1EpMTBCdyXyyoyOlodTpMVZZL7d35RoFpiIiIiFyL1WLB7kqN1J7pWQ+FgmtP9RUWVOSUZo+NlrYvaMOQdyVHXBgiIiIiF3Lp1EXkaRyz+82aNFxgGvfHgoqcUvzIIYg2lkntXcevCkxDRERE5Dp2nS+QthONRQiLDBWYxv2xoCKnJJPJMDvIJrX3NvrA3Mihf0REREQdqa+sxiGZYzKK2TEeAtMMDCyoyGmlTUmE0moGAFSovXDyyCnBiYiIiIic28GDp2FQ2J+f0psNSJ2cKDiR+2NBRU7Lx88HE2wlUnvn9WqBaYiIiIic364iq7Q9XV0BrVbTwdHUG1hQkVObM9zRZZ2pCkdFfqHANERERETOK/fiVWRrHc9LzRk/WGCagYMFFTm1MeNHIMBUAwCwyBXYe/ic4EREREREzmnXyevS9mBjGeLiB4kLM4CwoCKnplTIcZdXg9TeVaWD1WIWmIiIiIjI+ZgaDNhnajZVeig/5vcX/kmT05s1aYS0nacNQPYJ9lIRERERNXfi0ElUqTwBACqrCTOmcDKK/sKCipxeWHgQEk3FUnvnhSKBaYiIiIicz66bddJ2Csrg7e0pMM3AwoKKXMLswY5fCocUoaivqBCYhoiIiMh5lN28hZOaCKk9ZxQX8u1PLKjIJaSmJkJvNgAADAoNDqafFhuIiIiIyEnsOXoRVpn9Y32QqRqJY+MFJxpYWFCRS9BqVJiuqZLau0pssNlsAhMRERERiWc1mbC71jGSZ5avEQo5P+L3J/5pk8uYkxwrbWfrwnHrfLbANERERETiXcg4jQKtfXY/mc2KWZNGCk408LCgIpcRFxeFwaZyqb3rVI64MEREREROYGd2mbQ9xlKK4BD/Do6mvsCCilzKrDCltL3PHAhTfV0HRxMRERG5r7qSEhxWhkvt2XG+4sIMYCyoyKXMmJIAldW+sG+V2hMZB04KTkREREQkxv6DZ2FUqAEAnuYGpEzkcD8RWFCRS/H21CNV7uja3nnLIDANERERkRg2qxU7yxwjd9L0tVCrlB2cQX2FBRW5nLmJjq7t09oIFF3LEReGiIiISIBrp8/jui5Eas+dGCcwzcDGgopcTmJiHMJM9inUbTI5dh27LDgRERERUf/akZUvbQ8zlSImJkxgmoGNBRW5HJlMhtn+Zqm9u9EH5kajwERERERE/aehqhrpaNY7FaUVmIZYUJFLmjU1AQqrBQBQpvbBqSOnBCciIiIi6h8HD55Gg9JeROktBkyZkig40cDGgopckp+/DybYSqT2jus1AtMQERER9Z+dhVZpe5qqEjqtRmAaYkFFLmvu8CBp+4Q6DGV5hQLTEBEREfW9mxeuIlsbKrXnjh8kLgwBYEFFLmzM+BEINFYDAKwyBXYfPi84EREREVHf2nHqhrQdayxFXPwgcWEIAAsqcmFKhRyzfRqk9q4aD1gs5g7OICIiInJdjfUN2GcKkNpzQhUC01ATFlTk0mZPGgm5zT6OuEjji6yMc4ITEREREfWNI4dOo1alBwBoLEbMmDpacCICWFCRiwsKC8I4U5HU3pFdJjANERERUd/ZeatR2p4iK4WHl4fANNSEBRW5vLlDvKXtY4oQVJaWC0xDRERE1PtuX7+Fc5pmk1GMDheYhppjQUUuL3nSGPiZ7NOmm+VK7D3IYX9ERETkXnZmXJG2oxvLMCwxXmAaao4FFbk8pUqJu3SOdah2lqtgtVo7OIOIiIjIdRgbTdhj8JHacwItkMv5Md5Z8G+C3MKcFMe3NHkaP1w4dUFgGiIiIqLec/zwaVSr7M9LqawmpHEyCqfCgorcQtigSIxuLJDaO89xkV8iIiJyDztuOEbiTLIWw9vfV1wYaoUFFbmNubFe0vZhWQiqOTkFERERubjCG7dwRt1sMorEMIFpqC0sqMhtpEwaA29THQDAqFBh74EzghMRERER9cyOo5dgk9k/socbKzBq7HDBiejrWFCR21BrVJild3SJby9Xw2qxCExERERE1H1GgwG7Gv2k9lxORuGU+DdCbmXuJMe3NnnaAFw4niUwDREREVH3HTtwElUqTwD2ySjumjZGcCJqCwsqcivhUaEYYyqS2tsuFgtMQ0RERNR923MN0vZklMLH16uDo0kUFlTkdubHOdZpOKIMR2UBiyoiIiJyLbezryFLGy6154+NFJiGOqIUHaC57OxsrF27FpcvX4bZbEZkZCTmzZuHtLS0bl0vMzMTGzZsQE5ODgBg0KBBWLx4McaPH9+p8/fv34+//vWvAIBHH30U9913X7dyUP+akJoI36unUKnyhFmuxO5DZ/HgQ7NFxyIiIiLqtB3HrwGyaABAlLEcIxImCU5E7XGaHqqMjAwsX74cp0+fRkxMDMaOHYvCwkKsXLkS77//fpevt2XLFrz++uu4fPkyhg0bhlGjRuHatWt4/fXXsWXLljueX11djVWrVkEmk3Xn7ZBAKqUCc7zqpfbOaj0sJpPARERERESd11hXhz3mAKk9L1TGz6ROzCl6qGpra7Fy5UpYrVY899xzSElJAQBUVlbipZdewubNmzF+/HgkJCR06nr5+flYvXo1VCoVli9fjvj4eGn/iy++iNWrV2PcuHEIC2t/Hv/3338fBoMBU6dOxYEDB3r+JqlfzZkyEp/vKoZNJkeBxh9ZR09j7LQJomMRERER3dHh9JOoUQUBANQWE9KmcjIKZ+YUPVR79uxBfX09kpOTpWIKAHx9fbFkyRIAwKZNmzp9vS1btsBisWDOnDlSMQUA4eHhuP/++2GxWLB169Z2zz979iwOHDiABx98ECEhId14RyRaSGggkiwlUnvb5QqBaYiIiIg6x2azYUeeY2TNVGUZvLz0AhPRnThFQZWZmQkASE1NbfVaUlISVCoVsrKyYDQaO3W9kydPtnu9SZMmtfiZX2c0GvGvf/0LERERWLx4cad+HjmnecMDpe0MdQTKb+UJTENERER0Z7nnsnFB12wyiqQYgWmoM5yioMrNzQUAxMbGtnpNqVQiOjoaJpMJ+fn5d7xWXV0dSktLAdgnofi6gIAAeHl5oaSkBPX19a1eX7NmDYqKivCDH/wASqVTjIikbhqfPAKBJvtCvxa5ArsOnReciIiIiKhj20/mSNuDTeWIHzZIWBbqHOEVQ319Perq6gAA/v7+bR7j7++Pa9euobS0tM0iqbmmYsrDwwNarbbNYwICAlBTU4PS0lJER0dL+3NycrB582akpaVh5MiR3Xg3Ds8++2yrfWq1Gq+99hoAIDAwsNXr/ampWAwKChKao68tCLZh1Vej/XbUe+NH3t5QajRiQ7mAgXJ/UNfx3qD28N6g9vDe6Lz6igrsheNxk3uHeCE4OFhgor7nDveH8B4qg8GxYJmmnQ+6TfubH3un67V3rfauZ7Va8Y9//AN6vR7f+c537hycXML981Mht1kAACUaXxzcni44EREREVHbtm1JR71SBwDQWhqxYB6nSncFvdJD9eabb+LWrVtdOuepp55CXFxcb/z4Fmw2GwB0eWrJLVu24Nq1a3jyySfh5dXzVahXrFjR4eulpaVSVhGavgUoKSm5w5GuTS4HJthKcUxm/7Zn3flijJrk3u+5NwyU+4O6jvcGtYf3BrWH90bn2Gw2rM9tBOz1FKarK9FQX4uG+lqxwfqY6PtDJpN1OPN3Z/RKQVVSUtKp55uaa2xsBIAWw/IaGxuh17eexaStY9uj09nvwo56s75+vZKSEnz66acYMWJEtxcRJuc1f1Qojl20F6+Z2kiUXMtB0JBBYkMRERERNXP91Dlc0Tk+2M+fMERgGuqKXimomp4L6g69Xg+9Xo/6+nqUl5e3WVCVl5cD6NxzR03H1NXVwWAwtFmElZWVtTj2/PnzaGxsRHV1NV5++eUWxzZVy7t27cLp06cxfPhwPPLII114hyTamLFDEXL2BIpU3rDK5Nh59BIeZUFFRERETmT7mduAejAAYKipDEOGDBeciDpL+KQUABATE4OLFy/i+vXriIyMbPGa2WxGbm4uVCoVwsPD27mCg4eHBwIDA1FaWoqcnBwMH97yZiwrK0NNTQ0CAwNbFW95eXnIy2t7au3i4mIUFxe3WfCRc1PI5ZgTaMEHVfb2TqM/vtHQANVXvZlEREREItWVlmG/3PE5d14MP6O4EqcoqJKSknDx4kUcPXoU06dPb/HayZMnYTKZMG7cOKjV6k5fb8eOHTh69GirgurIkSPSMU3S0tLaHeq3Zs0afP7553j00Udx3333df5NkVOZPTURn2zMgVmuRLnaGxnpJzBl3jTRsYiIiIiwL/00DMoIAIDeYsC0yYmCE1FXCJ/lDwBmzZoFnU6HEydO4NixY9L+qqoqfPDBBwCARYsWtTpv2bJlWLZsmTQksMmCBQsgl8uxc+dOXL58WdpfUFCAdevWQS6XY8GCBX30bsgZ+fl6YrK8TGpvvdW5RaKJiIiI+pLVYsbWCkenwSxdNbQalcBE1FVO0UPl6emJJ598Em+99RZWrFiBkSNHwsvLC1lZWairq8Pdd9+NxMTWlXrTRBhms7nF/vDwcCxZsgSrVq3C8uXLMXr0aCgUCpw9exZGoxFLly7t1PBBci93j41C+ml7IZWli0DuhcuIHhkvOBURERENZBcyzuKW1rEG0/zUYQLTUHc4RUEFAKmpqXj55Zexdu1aXLlyBWazGREREZg3bx5mzpzZ5estWrQIoaGh2LhxIy5evAgAiI2NxeLFi5GcnNzb8ckFjBg5GDEZR3BTbV9AeltmDn7IgoqIiIgE2nKxFNBEAwBGm4oQGcXJKFyNzCZyMaQBrKCggOtQCbB1+zH8vdQHAKA3G/Cf+2Oh8/YWnMr5DNT7g+6M9wa1h/cGtYf3RvvK8wvxxO4yWOQKAMAvBzdi0uQxglP1L9H3R2+sQ+UUz1AR9ZcZ08ZAZ7GvQ1av1GL//tNiAxEREdGAtfPgeamY8jfVYMLEBMGJqDtYUNGAotdrMVPlmMRkW7EMVqtVYCIiIiIaiMxGE7bXeUnteT4NUCoVAhNRd7GgogHn7olx0vYNbRCyT18QmIaIiIgGouOHT6FMbX/sQGG1YM5U9k65KhZUNOBED4lCgrFIam/JKhSYhoiIiAairddrpe1UaxECgvwFpqGeYEFFA9KCwXpp+7A8FBUlZR0cTURERNR7bl/PxRmNYwmfBQmhAtNQT7GgogFp4uSx8DPVAADMciV2HTwnOBERERENFNuOXZG2oxrLMHIcp0p3ZSyoaEBSqVWYq3d0tW+v1MJstghMRERERAOBob4Be4yO4X13B1shl/MjuSvj3x4NWHOnjoDcZi+iStQ+OHnsrOBERERE5O7S00+hTqkDAGgtjUibNlZsIOoxFlQ0YAWGh2KiyTEhxdbsCoFpiIiIyN1ZrVZszXeMiEmTl8DDx6uDM8gVsKCiAW3BqGBp+5Q6FPnXcwWmISIiInd25cwlXNcESe35E4YITEO9hQUVDWiJyaMQYbT3TNlkcmw7ki04EREREbmrLWfypO2RjUUYPGywwDTUW1hQ0YAml8sxP9gqtXeZA2Goqe3gDCIiIqKuqygoxkFFmNS+O9ZTYBrqTSyoaMCbNWMctJZGAECdUod9+zIFJyIiIiJ3s+PAWZjlSgCAn6kWk6aOEZyIegsLKhrwPDz1mKlyTEixpVgOq4VTqBMREVHvMDU2Yludt9Se51MPlVIpMBH1JhZURAAWpA6Vtm9qg3A+g1OoExERUe84uv8EytX2gkppNWPe9NGCE1FvYkFFBCB6cARGm4qk9uaLpQLTEBERkbuw2WzYktsotSejBP4BvuICUa9jQUX0lUXDHKuWH1OHozjnlsA0RERE5A5ysi7igi5cai8cHyMwDfUFFlREXxk/cRSCjVUAAKtMgW2HLgpORERERK5u00nHF7RDjGUYNoJTpbsbFlREX1Eq5Jgf6JiMYqcpAI21dQITERERkSurLipButwxVfrCGA1kMpnARNQXWFARNTMnbQzUFhMAoFrlgYP7TghORERERK5qV/ppGBVqAIC3uR5TOVW6W2JBRdSMt5cHpivLpPbmQnAKdSIiIuoyc2MjttZ4Se05nnXQqFUCE1FfYUFF9DULm02hfk0XgisnOIU6ERERdU3mgRMo1vgCAOQ2C+ZPTxAbiPoMCyqir4mNjcBIU4nU3nS+WGAaIiIicjU2mw2bcxqkdoqtFMFBfgITUV9iQUXUhoVDfaTtw6oIlOfeFpiGiIiIXMnt89k4o4uU2gvHRnZwNLk6FlREbUhJTYC/qRYAYJYrsf3gecGJiIiIyFVsPnFD2o4xlSMhIVZgGuprLKiI2qBSyDHf3yi1tzf6wVTPKdSJiIioY3UlJdgrcyzkuyCKU6W7OxZURO2YO300lFYzAKBC7Y3DnEKdiIiI7mD3vlMwKDUAAA+zATOmjhaciPoaCyqidvj5emKqolxqb8y3wma1CkxEREREzszcaMTmak+pPduzFjoNp0p3dyyoiDpwz0THmOcrujBkcwp1IiIiakfmgRMo1PoDAOQ2KxZOHSk4EfUHFlREHYiLi2wxhfrGc5xCnYiIiFqz2WzY2GKq9BKEhPgLTET9hQUV0R0siveWtg+rI1DCKdSJiIjoa3KyspGli5Dai8ZGdHA0uRMWVER3kJKSiCBTNQDAKlNg68GLghMRERGRs9l0KlfajjWWYlTCEIFpqD+xoCK6A6VCjgUBJqm93RgAQy2nUCciIiK7yqJS7JeHSe17YjhV+kDCgoqoE+bMGAOtxb4uVa1Kj337TgpORERERM5ie/pZmOT22fx8TbWYOnWs2EDUr1hQEXWCl7cn0pSlUntTsRxWTqFOREQ04Bkbjdha6yW153vVQa3mVOkDCQsqok66Z9JQafuWJgBnM7IEpiEiIiJncDj9JCrU9oJKaTVj/nQu5DvQsKAi6qTIwVFIMuZL7Q0Xyzs4moiIiNydzWbDxltGqT0dRfAL8hOYiERgQUXUBfcMd6wnkakOQ94NTqFOREQ0UF3KuoyrmmCpvWh8jMA0JAoLKqIuGDsxEZGNZVJ705HLAtMQERGRSBtPOb5YHdVYiCEj4wSmIVFYUBF1gVyhwKJgx2QUe8wBqKmoFpiIiIiIRCi+XYAjilCpfU+sh8A0JBILKqIumjkzGZ7mBgCAQaHBrr2ZghMRERFRf9uSfh5WmQIAEGyswoSp4wQnIlFYUBF1kdZDhzm6Sqm9pVIHs9HY/glERETkVhqqa7DT5HiuekGACUqlUmAiEokFFVE3LJieCLnNAgAo1vji6P4TghMRERFRf9mzJxO1Sj0AQGsxYnYae6cGMhZURN0QHBqIybZiqb3hphE2m01gIiIiIuoPZqMJG8s1UnuWuhxe3nx+aiBjQUXUTfclD5K2s3WhuHTirLgwRERE1C+OHziBAo19rSm5zYrF00cKTkSisaAi6qahIwZjpMnRS7U+q0hgGiIiIuprNpsNG240SO0UazFCw4M7OIMGAhZURD1w7zDHaujH1JEovHpDYBoiIiLqS1dOnsMFXbjUvnd8lMA05CxYUBH1wISUUQgzVgIArDI5Nh7mQr9ERETuav2ZAml7mKkEI0YNEZiGnAULKqIeUMjluCdcJrV32UJRU1IqMBERERH1heIbN3FYHSG17x3qLTANORMWVEQ9dNeMcY6FfpUa7Nh7SnAiIiIi6m2bDl6SFvINMVUhJTVRcCJyFiyoiHpIp1Vjnlet1N5c6w1TQ0MHZxAREZErqSsrw05riNReGGKDUsGP0WTHO4GoFyxMGw2l1QwAKNP44NCeDMGJiIiIqLfs3HMS9UotAEBvNmAOF/KlZlhQEfWCAH8fTFOUSe31BYDVYhaYiIiIiHqDucGATTVeUnuuRw30Ok0HZ9BAw4KKqJfcOylO2r6uC8H5wycFpiEiIqLecGRvBko0vgAAhdWCRWl8dopaYkFF1EsGD47AGHOJ1F6fXSEwDREREfWU1WLB+nyr1J4iL0NQoK+4QOSUWFAR9aJ7Ex2rpZ/QRuH2uUsC0xAREVFPXDp2Eld0oVJ7cWqswDTkrFhQEfWipHHxiDLZe6ZsMjk2HL8hOBERERF115cXyqXtUeZSDB0SKTANOSsWVES9SCaTYXGM40HVvYpIVNzOF5iIiIiIuuP2hWxkaKOk9r0jAwWmIWfGgoqol82YOhq+5joAgFGhwpb9WYITERERUVetz7gBm8z+UTnCVIkJycMEJyJnxYKKqJdpVEos8m+U2ltMQWioqhKYiIiIiLqiIq8AexURUvu+aCXkMpnAROTMWFAR9YH5aWOhtdiLqlqVHrt2ZwpORERERJ21aX8WTHIVAMDXVIu0aWMEJyJnxoKKqA94eekxR1sptTdU6mE2GsUFIiIiok6pr67GNpPjealF/o1Qq1QCE5GzY0FF1EcWTx8FhdUCACjW+OLQvhOCExEREdGd7NydiVqlHgCgtTRi/syxYgOR02NBRdRHgkMDMQ3FUnvdLQusVmsHZxAREZFIJqMRGyr1UnuepgJeXh4CE5ErUIoO0Fx2djbWrl2Ly5cvw2w2IzIyEvPmzUNaWlq3rpeZmYkNGzYgJycHADBo0CAsXrwY48eP7/C8rKwsbN26FVeuXEFdXR28vLwwaNAgzJkzB8nJyd3KQgPTfamx2JfRAAC4oQ3C2YwsjE3lOGwiIiJndHB/JkrVfgAAhdWCRdNHCU5ErsBpCqqMjAysWLECNpsNI0aMgJeXF86dO4eVK1fi5s2beOyxx7p0vS1btuC9996DQqFAYmIilEolzp49i9dffx2PP/44FixY0OZ5H374IdavXw+lUolhw4bBx8cHFRUVuHDhAvz8/FhQUZcMHhqDcQf34pQ6DACw7kI5xqYKDkVEREStWK1WfHnLAny1nOQ0FCI4jAUV3ZlTFFS1tbVYuXIlrFYrnnvuOaSkpAAAKisr8dJLL2Hz5s0YP348EhISOnW9/Px8rF69GiqVCsuXL0d8fLy0/8UXX8Tq1asxbtw4hIWFtThvx44dWL9+PYYMGYLnnnsOgYGOBxIbGxtRVFTUS++YBpL7RwXg1BX79mlNGK5duo4hw2PFhiIiIqIWTh8/hxyN47PffRP5/2rqHKd4hmrPnj2or69HcnKyVEwBgK+vL5YsWQIA2LRpU6evt2XLFlgsFsyZM0cqpgAgPDwc999/PywWC7Zu3drinLq6Onz44YfQ6XT4xS9+0aKYAgCNRoPo6OjuvD0a4BKTEzCk0fEs1ZcZOeLCEBERUZvWXSiXtpOM+Rg8bLDANORKnKKgysy0r9GTmtp6LFRSUhJUKhWysrJg7OS00ydPnmz3epMmTWrxM5scOnQIDQ0NmDJlCvz8/LqUn6gjcrkc90erpfZBeQiKbuULTERERETNXT1/FWfVoVL7/pGBHRxN1JJTFFS5ubkAgNjY1l2rSqUS0dHRMJlMyM+/84fQuro6lJaWArBPQvF1AQEB8PLyQklJCerr66X9WVlZAIDRo0ejsrISmzZtwj//+U+sXr0aGRkZnJ2NemTSjGSEGCsBAFaZAhvSL4gNRERERJJ1x3Ok7SGNxUiY0LnHTIgAJ3iGqr6+HnV1dQAAf3//No/x9/fHtWvXUFpa2maR1FxTMeXh4QGtVtvmMQEBAaipqUFpaak0jO/27dsAgJKSEvz9739vUWxt3LgRgwcPxvPPP99uxq979tlnW+1Tq9V47bXXAKDVkML+plTa/+qDgoKE5hhIHo5S4M9fPYa30xKMJ602+IUEiw3VDt4f1B7eG9Qe3hvUHme/N3IvX8NhZbjU/nZCAEJCQgQmGlic/f7oDOE9VAaDQdrWaDRtHtO0v/mxd7pee9dq73q1tbUAgI8++gihoaH4/e9/j/fffx+/+93vMHjwYNy4cQNvvvkmbDbbHTMQteW+xWnwMtsL9UaFGp9+mS44EREREX249TisMvtH4lBTFWbPmyY4EbmaXumhevPNN3Hr1q0unfPUU08hLi6uN358C00Fj0wm69J5TUP61Go1XnjhBXh7ewMA4uPj8cILL+Cpp57ClStXkJWVhdGjR9/xeitWrOjw9dLSUqHFWdO3ACUlJcIyDEQLvWvxSb19wcAvqz2x8No1aL+615wJ7w9qD+8Nag/vDWqPM98bFflF2G4OBhT29uJQGyoqyjs+iXqV6PtDJpO1mvm7q3qloCopKenU803NNTY2AkCLYXmNjY3Q6/WdOrY9Op0OQMe9WW1dT6fToaamBuPHj5eKqSY+Pj5ISkrCkSNHcOHChU4VVERtWXDXOHy5/gYMCjWqVR7Yses4Fj8wS3QsIiKiAWnjvjMwKuyPf/iY6jArLUlwInJFvVJQNT0X1B16vR56vR719fUoLy9vs6AqL7d/U9CZ546ajqmrq4PBYGizCCsrK2t1vaCgIBQXF7c7frNpf1VV1R0zELXHx8cLc7QV2Giyj81eX+mJ+Q0GqHV3/rKAiIiIek9dWTm2moKkT8P3+Bug1ao7PomoDcKfoQKAmJgYAMD169dbvWY2m5GbmwuVSoXw8PBWr3+dh4eHVCjl5OS0er2srAw1NTUIDAxsUbwNHmxfa6DpWaqva9rfmV4yoo7cO3M0lFYzAKBU44P0XUcFJyIiIhp4tu7ORL3SPrJJbzHg7lnsnaLucYqCKinJfgMfPdr6g+XJkydhMpmQkJAAtbpz3xp0dL0jR460OKZJcnIyAODChQutpki3Wq24ePEiAEfhRdRdQUF+mKEok9pri1WwmDq3xhoRERH1nKGmBhvrHOuO3u1RA08PncBE5MqcoqCaNWsWdDodTpw4gWPHjkn7q6qq8MEHHwAAFi1a1Oq8ZcuWYdmyZdKQwCYLFiyAXC7Hzp07cfnyZWl/QUEB1q1bB7lcjgULFrQ4Z+TIkYiPj0deXh7Wrl3b4rXPPvsMBQUF8PHxwcSJE3v8fokemDYcMpu9cM/TBuDo3gzBiYiIiAaOPbuOo1LtCQBQW024566xYgORSxO+DhUAeHp64sknn8Rbb72FFStWYOTIkfDy8kJWVhbq6upw9913IzExsdV5TRNhmM3mFvvDw8OxZMkSrFq1CsuXL8fo0aOhUChw9uxZGI1GLF26tM3hg0899RRefPFFrFmzBocOHUJkZCRu376NvLw8qNVq/OQnP+GQP+oVkVEhSMUlHIH9Waovci2YZDFDrnCKf5JERERuy2ww4MsKHfDVCjuz1BXw82v9OZOos5zm01tqaipefvllrF27FleuXIHZbEZERATmzZuHmTNndvl6ixYtQmhoKDZu3CgN14uNjcXixYul4X1fFxoaijfeeANr1qzBqVOncOLECXh6emLKlCl44IEHEBUV1aP3SNTcQymxOJJhX9T6mi4EZw9mYuyMFMGpiIiI3NuB3cdQpLF/oSm3WXBf2ijBicjVOU1BBQDDhw/HCy+80Onj16xZ0+HrycnJ7RZP7fH19cUPf/jDLp1D1B1xQ6Mw9nA6TiuDAQCfX6nBmOm2Lq+hRkRERJ1jMZmwtkAOfPW41HR5KUJDWFBRzzjFM1REA9WDSRHSdpYuEpePnxYXhoiIyM2d2J+BXJ1jiZwHpsQLTEPuggUVkUCJCbGIN5VK7c/PFgtMQ0RE5L6sFgu+uNEotSdaixATEyYwEbkLFlREAslkMjw0KkBqZ+hikHv2gsBERERE7unCkUxk6x2Tkj00cZC4MORWWFARCTYheTiiTBVS+4vjuQLTEBERuR+bzYbPL1ZK7QRzCYYNixEXiNwKCyoiweQyGR4copfa6ZooFF6+JjARERGRe7meeQan9NFS+6GxoQLTkLthQUXkBKZNTkSwqRoAYJUpsO7w5TucQURERJ215nSRtD3EVI6xo+MEpiF3w4KKyAkoFXI8EOn457hLEYWym7cEJiIiInIPN7Mu4ajOMbzvGyN8uEQJ9SoWVERO4q4ZY+BvqgUAmOVKrEvn5BREREQ99VmzZ5OjTRVImThSYBpyRyyoiJyERqXCfSFmqb3dFoaKAk6jTkRE1F15l6/hkDpSan9jiA5y9k5RL2NBReRE5s0cBx9THQDAqFBj/d4zghMRERG5rs+PXINVZv+4G26sxOTJowUnInfEgorIiWi1Gtzr1yC1t5qDUV1W0cEZRERE1JbCm3nYL4+Q2g9GyaFU8KMv9T7eVURO5u5Z4+FprgcAGBQabNx9SnAiIiIi17P2wEVY5AoAQLCxCjPSkgQnInfFgorIyeg9dbjHs1pqb27wQ211rcBERERErqW0sBi7bY61ph4IsUClVApMRO6MBRWRE1o4Kwl6swEAUKfUYevuTMGJiIiIXMe6vedgltsLKH9jDWbNShaciNwZCyoiJ+Tl640F2jKpvaHKEw119QITERERuYaKknLsMAdJ7fsD6qHWqAUmInfHgorISS2+axy0lkYAQLXKA9t2HheciIiIyPmt330KRrkKAOBjqsXc2RMEJyJ3x4KKyEn5BPljnqpUan9Z6QFDPXupiIiI2lNVWo6txkCpfa9PLbR6vcBENBCwoCJyYvfNGgu1xQQAqFR5YteOY4ITEREROa+NO0/CoNAAADzN9bh77kTBiWggYEFF5MT8gwMwR+3opVpXroeRz1IRERG1UlNWhs3GAKm92LsWeg/2TlHfY0FF5OTumzkaSqsZAFCq8cGuHUcFJyIiInI+m3Zkol6pAwDoLQYsmD1ecCIaKFhQETm54JAAzFY6eqm+qPCAsa5OYCIiIiLnUlNahg1Gx8x+93jVwMvLQ2AiGkhYUBG5gAdnjWEvFRERUTs27Tzh6J0yG3APe6eoH7GgInIBwcF+mK1yrEv1eaUne6mIiIgA1JSUYoMxWGrf410LLy8+O0X9hwUVkYt48C7Hs1Rlah/s2n5EcCIiIiLxNu7MbNk7NSdJcCIaaFhQEbmI1r1UXjDWspeKiIgGrpqSEmw0Neud8qmDlyd7p6h/saAiciEPzmrWS6VhLxUREQ1sG3ec/NqzU+MEJ6KBiAUVkQsJDvLDbHW51P68ygvGmlqBiYiIiMSoKS7BRrOjd2qxL3unSAwWVEQu5sG7Elv0Uu3cwV4qIiIaeFo8O2UxYNFsPjtFYrCgInIxwUF+mK1p1ktV7QNjTY3ARERERP2rpqgYG80hUnuxTx28PHQCE9FAxoKKyAU91GzGv3K1N3Zu57pUREQ0cGzY5Xh2ysPM3ikSiwUVkQsKCvTFnOa9VDU+aGQvFRERDQA1RUXYZA6V2vf41bN3ioRiQUXkoprP+Feu9sbObeylIiIi97dh50nUK7UAvuqdmsXeKRKLBRWRiwoK8MUcTYXU/qLWB4Zq9lIREZH7qiooxCZLmNS+x68BXh5agYmIWFARubSHZo2Gqlkv1TY+S0VERG7sy92nW/RO3TOL606ReCyoiFxYYIAP5msdz1J9UeeH+qoqgYmIiIj6Rnl+ITZbw6X2/f4N8GTvFDkBFlRELu7BWWOhsRgBANUqT2zaflxwIiIiot73xe6zaFSoAQDe5nos5Mx+5CRYUBG5OD9/byzUO56l+rIhADUVleICERER9bLiW/nYBkfv1IOBBuh1GoGJiBxYUBG5gfvnJEFvNgAA6pQ6rN+RKTgRERFR7/ls3wWY5UoAgL+pBnfPThaciMiBBRWRG/D28cJiz0qpvbExEFWlFe2fQERE5CIKbuZhNxwz+30j2ASNRi0wEVFLLKiI3MQ9c5Phaa4HABgUGqzddUpwIiIiop77NP0SLHIFACDYWIXZsyYITkTUEgsqIjfh6eWJ+70cM/xtMQairLhMYCIiIqKeuXXjFvbLHL1T3wwzQ61RCUxE1BoLKiI3snDuBPiaagEARoUaX+w6LTYQERFRD3ySfhlWmf3janhjOdLumig4EVFrLKiI3IjO0xMP+tZJ7e2WEBTn5gtMRERE1D3Xz1/FQWWE1H4kSg6lmr1T5HxYUBG5mXnzUxBgrAYAmOVKfLr3nOBEREREXffRsRxpO6axDFPv4rNT5JxYUBG5GY1Wi4fDLFJ7jyIC+VeuC0xERETUNdknzuK4JlJqPxqvh0KhEJiIqH0sqIjc0KzZExBitE9QYZUp8NGBK4ITERERdY7NZsMHp4ul9hBjKSZOHiMwEVHHWFARuSGVUolHBiul9gFNDK6d5tA/IiJyfmcOHMdZnaN3asnoQMjl/MhKzot3J5Gbmj5tHKJMjsV9PzhRAJvNJjARERFRx6xmM1ZdbpDaCaYSjEsaJjAR0Z2xoCJyU0qFHN8Z4S21T+qicO7QCYGJiIiIOnZ412Fc04VI7e+kRkEmkwlMRHRnLKiI3NjEiSMxzFQqtVdfqoXVbBaYiIiIqG3mBgM+zHNMPJFiLcLw4YPEBSLqJBZURG5MJpNhaXK41M7WhSFj71GBiYiIiNq2e9sh5GsDAABymxVLpscLTkTUOSyoiNxcQkIsxlscsyV9kGuDubFRYCIiIqKWDNU1+KTSS2qnKUoRHRMmMBFR57GgIhoAvjN1CGQ2KwDgljYI+7YfEpyIiIjIYfOWwyhX25/7VVnN+NasBMGJiDqPBRXRADA4NgLT5I5nqT4u1cNYUyswERERkV1NUTHWGoKk9t26CgQH+wtMRNQ1LKiIBohHZ46E0mqfkKJU44utWw8LTkRERASs25GJWpUeAKC3GPDQ3CTBiYi6hgUV0QARFhaIuRrHulSf1wegrqS0gzOIiIj6VllOLjbaIqT2fb4N8PH2EJiIqOtYUBENIA/PHQutxQgAqFZ5YP3244ITERHRQPbp3nMwKtQAAB9zHe6ZM15wIqKuY0FFNID4+XrhHm/Hs1PrLeGouHVbYCIiIhqo8s9fwi5ljNR+OMwGvU4tMBFR97CgIhpg7pubBC9zAwDAoNTgk11ZghMREdFAY7PZsOrITVjk9oV8Q0zVmDtznOBURN3DgopogPHUa/FQiElq71BF4/bFKwITERHRQJOdcQZHdI7eqUfjtFCrFAITEXUfCyqiAWjBrCSEmKoAAFaZAqsO54gNREREA4bVYsF/sxyTJA0xlWH6lESBiYh6hgUV0QCkVinxnSEaqX1MG4Xzx88KTERERAPF0X3HcUkXJrUfHxcEuUwmMBFRz7CgIhqgpkwZjTijY9r0985WwGq1CkxERETuzmQ0YtVNm9RONhdi9Jh4gYmIeo4FFdEAJZfL8d2xgVL7sjYEh/efEJiIiIjc3fYdx1Cg8QMAyG1WPDYlVnAiop5Tig7QXHZ2NtauXYvLly/DbDYjMjIS8+bNQ1paWreul5mZiQ0bNiAnJwcAMGjQICxevBjjx7e9xoHVasXOnTuxf/9+3L59GyaTCX5+fkhISMD999+PsLCwNs8jclUJ44Zjwpk9OK4KBwCszrFiYqPpDmcRERF1XV1NHT4t1QMqe3sWChAdN0tsKKJe4DQ9VBkZGVi+fDlOnz6NmJgYjB07FoWFhVi5ciXef//9Ll9vy5YteP3113H58mUMGzYMo0aNwrVr1/D6669jy5YtrY632Wz405/+hHfffRe3bt3C8OHDMWHCBCgUCuzbtw/PP/88rl271htvlcipPDZ1COQ2CwCgUO2LbTuOCU5ERETuaO2246hWeQAAtJZGfGs2J6Ig9+AUPVS1tbVYuXIlrFYrnnvuOaSkpAAAKisr8dJLL2Hz5s0YP348EhISOnW9/Px8rF69GiqVCsuXL0d8fLy0/8UXX8Tq1asxbty4Fj1OmZmZOHHiBIKDg/G73/0Ovr6+AOy9VqtXr8bmzZuxatUqvPzyy7375okEi4qLwZzDO7FdEQUAWFOmx0NlFfAJ8BOcjIiI3EVhbh42NPgDX82Mfq+mFAGhY8SGIuolTtFDtWfPHtTX1yM5OVkqpgDA19cXS5YsAQBs2rSp09fbsmULLBYL5syZIxVTABAeHo77778fFosFW7dubXHOhQsXAACzZ8+WiinA/pzJgw8+CADsoSK39cjcsdBaGgEANSo9/vvpbsGJiIjInfzti4MwKtQAAF9TLe5bkHKHM4hch1MUVJmZmQCA1NTUVq8lJSVBpVIhKysLRqOxU9c7efJku9ebNGlSi5/ZRKVStXs92VdTeXp6enbq5xO5Gv/QINzn4VgTZJ3BH3lXrgtMRERE7uLiiTPYhVCp/UiwEXovfqYi9+EUBVVubi4AIDa29UwvSqUS0dHRMJlMyM/Pv+O16urqUFpqnwp60KBBrV4PCAiAl5cXSkpKUF9fL+0fPXo0AGD37t2orKyU9lutVnz22WcAgBkzZnT6PRG5mvvuToGfqRYAYJKr8Pf1RwUnIiIiV2ez2fDXXRdhk9k/ckYayzFn7kTBqYh6l/BnqOrr61FXVwcA8Pf3b/MYf39/XLt2DaWlpW0WSc01FVMeHh7QarVtHhMQEICamhqUlpYiOjoaADBq1CgsWrQImzZtwk9+8hOMGDECWq0WN27cQHl5ORYsWICHH3640+/r2WefbbVPrVbjtddeAwAEBga2er0/KZX2v/qgoCChOci5fDfuKlbctG/vUUTgG9nXMXoqh2WQA393UHt4b1Bb0jftxklthNR+MjmUsyZTC+7wu0N4D5XBYJC2NRpNm8c07W9+7J2u1961Orre0qVLsXTpUpjNZpw+fRpHjx5FUVERwsPDMXLkSCgUijv+fCJXdt/iNMSY7EP/bDI53jl4E1azWXAqIiJyRaaGBvz1XLXUTrSWIW0me6fI/fRKD9Wbb76JW7dudemcp556CnFxcb3x41uw2eyrbzc999RZJpMJf/nLX3Ds2DE88MADSEtLg5eXF65du4b//ve/+NOf/oTvfe97mD9/fqeut2LFig5fLy0tlbKK0PQtQElJibAM5JweHxOAly9YAQAXtGHY8PF6TJk/XXAqchb83UHt4b1BX7f5i124pY0EAMhsVnx3Uow0koioiejfHTKZrMe9pr1SUJWUlHTq+abmGhvtM4o1H5bX2NgIvV7fqWPbo9PpAHTcm9XW9datW4cjR460GtqXkJCAX/3qV/jpT3+Kjz76CFOnTuXkFOTWksbFY8KFgzgO+7DU9/NVSK6phYYPEBMRUSfVlJTg4xo/aRHfuZoKDIkbKTYUUR/plYKq6bmg7tDr9dDr9aivr0d5eXmbBVV5eTmAzj131HRMXV0dDAZDm0VYWVlZq+sdOHAAQNszAwYGBiI+Ph5ZWVm4du0axozhugnk3p5ZnIzHvrwGi1yBIo0fNm0+hAcfmSc6FhERuYg1W06gRjUYAKC1GPHkN2YA4BByck/Cn6ECgJiYGADA9eutp2k2m83Izc2FSqVCeHj4Ha/l4eEhFUo5OTmtXi8rK0NNTQ0CAwNbFG9NRVZbBR3g6Pmqra29YwYiVxcbG4V7PB3j3j9rDEHFrdsCExERkavIu3gZW+RRUvvRUDOCg7hYPLkvpyiokpKSAABHj7aepvnkyZMwmUxISEiAWq3u8fWOHDnS4pgmTYv5trV4r9VqxY0bNwC49gwkRF3xo2/eBU9zAwCgQanFh7uyBCciIiJnZ7PZ8N6hGzDL7YOgAk01WPLQTMGpiPqWUxRUs2bNgk6nw4kTJ3Ds2DFpf1VVFT744AMAwKJFi1qdt2zZMixbtkwaEthkwYIFkMvl2LlzJy5fviztLygowLp16yCXy7FgwYIW50yYMAEAsGbNmhbPg1mtVnz00UcoKSlBUFAQhgwZ0vM3TOQCfHw88c1Qk9TerYrB9ZMsqoiIqH1n0jOQoYuR2o8P1UKr6dwX4kSuSvg6VADg6emJJ598Em+99RZWrFiBkSNHwsvLC1lZWairq8Pdd9+NxMTEVuc1FT7mr03rHB4ejiVLlmDVqlVYvnw5Ro8eDYVCgbNnz8JoNGLp0qWthg8+9NBDOHPmDPLz8/Hzn/8c8fHx8PT0RE5ODoqKiqBWq/Hkk09y6nQaUO6enYxtqzOQp/KFVSbHfzLz8cqYEZArnOJXBxERORFzowH/udII2J+SwDBzKaZOniI2FFE/cJpPRampqXj55Zexdu1aXLlyBWazGREREZg3bx5mzux6V/GiRYsQGhqKjRs34uLFiwCA2NhYLF68GMnJya2O9/LywquvvoqNGzfi+PHjuHr1KsxmM/z8/DBjxgzce++9iIyM7PH7JHIlKoUc3030xe8u2dtZ+khk7DqC1HnTxAYjIiKns2vzQdzUOT4rfX9STJeXsSFyRTKbyMWQBrCCggKuQ0VOq/n9YbPZ8P9WH8RphX1fmKEc7zwyBmoPD5ERSRD+7qD28N4Y2OpKSvHkphxUqe1LbMyQF+PZb9nXMOS9QR0RfX/0xjpUTvEMFRE5L5lMhu/NGAq5zb7Yb4HWH5s3HRScioiInMlnWzOkYkptMeE7c0cLTkTUf1hQEdEdxcSEYq66TGqvMYSg4nbXFvMmIiL3lHfpCjbKoqX2A761CArwFReIqJ+xoCKiTnl0/nh4mA0AgHqlFqt3csY/IqKBzmaz4d1DOdI06QGmGtw3r/Wz6kTujAUVEXWKj7ce3wo1Su3d6hhkZ54TmIiIiEQ7kX4cmVrHIr7fHaqBTqMSmIio/7GgIqJOu3vWeEQbHeu+/et0GSwWcwdnEBGRuzIaGvHuNYvUHmUqwdTJrZe5IXJ3LKiIqNOUSgV+MMZfal/RhmDPzgyBiYiISJT1W46gQOMHAJDbrPjBFE6TTgMTCyoi6pLRScMxxeyYkGJ1oQa1NbUCExERUX8rLSzF57V+Unu+vBCDh0Z3cAaR+2JBRURd9vjskVBb7M9TVak88Mnm44ITERFRf3p/51kYFBoAgLepDo/ezYkoaOBiQUVEXRYcEYqHdKVSe7M5GLlXcsQFIiKifnPu5AWky8Ol9rcD6+Hl5y0wEZFYLKiIqFvuWzgJIcZKAIBVpsC/DtyA1WoVG4qIiPqU2WjCv047JieKNZRg9rxJAhMRiceCioi6RaPX4XtDHFPjntWE4cieYwITERFRX9ux9RByNIFS+4fjg6BUKQUmIhKPBRURdVvKtCSMNRVK7f/kymGoqhaYiIiI+kpVYTE+rHAM7Uuz5mPE2OECExE5BxZURNRtMpkMP5gZD4XVvg5JqcYHn286LDgVERH1hQ+3ZqJWpQcA6CyNWHr3OMGJiJwDCyoi6pHImHAs0jnG06+zRiLvwmWBiYiIqLdlZ5zCDlWM1H44oAEBgX4dnEE0cLCgIqIee2ThRPib7GtRmeVK/PPQTVgtZsGpiIioN5gbDfj72WrYZPaPjVGmCiyaO0FwKiLnwYKKiHpMr9Pg+8M0Uvu0PgqHth8UmIiIiHrLtg3puK4Lkdo/Sg6GWqUQmIjIubCgIqJeMSU1AWMtJVL7P0V61JeVCUxEREQ9VZ57Gx/WB0ntNHkxEhOGCExE5HxYUBFRr5DJZPjRnOFQWe1D/crV3vh4E6dRJyJyVTabDe/tzEK9UgcA8DAb8DgnoiBqhQUVEfWa8LAgPOBdI7U3KWJw42SWwERERNRdWekZ2K8dLLWXRJjh5+slMBGRc2JBRUS96oG7kxFqqgIAWGUK/C2zBBaTUXAqIiLqCmNdHf5xxSS1h5jKMHdmksBERM6LBRUR9SqtWoUfjHFMpZutD8fuTekCExERUVet33AAt3WBAACZzYonp0RDqeDHRqK28F8GEfW65HHxmATHBBWrqvxQlV8gMBEREXVW0eVrWGMKl9rz1GUYOjRKYCIi58aCioj6xPfmJUJrsQ/1q1F5YPW2U4ITERHRndisVvx7/1UYFWoAgI+5DksWJgtOReTcWFARUZ8IDvTFI8EGqb1TE4uLh08ITERERHeSsfMQMvQxUvu7sSp4eegEJiJyfiyoiKjPLJqTjGhThdT+68UGGBsaBCYiIqL21FdU4h/5aqk9ylKKtKmJAhMRuQYWVETUZ1QKOX6cGia1b2mDsG7DQYGJiIioPR9sPIYytQ8AQGk1439mDoVMJhOcisj5saAioj41YvggzFcUSe3PjKG4fTVHXCAiImol+/QFbJE7Jp54yKsS0VEhAhMRuQ4WVETU576zaCL8TLUAAJNchb+l58BqtQpORUREAGAymvDXk+WwyewfCyOMFXhwQYrgVESugwUVEfU5T08dfhjrGDZyThOK3bsyBCYiIqIm67ccwU1NoNT+8Vg/qNUqgYmIXAsLKiLqF5OmJWGiMU9qv1egQUVpucBERESUfzMPn9Y4FmOfa72NhHHDBSYicj0sqIioX8hkMvxw7ihoLY0AgFqlDv/ZzLWpiIhEsVqt+NvuyzDK7b1RvsYaLF3ENaeIuooFFRH1m6CIUHzHt1JqpysjkHmAa1MREYmwb8dhnNU4ZmJ9ItoCLz9fcYGIXBQLKiLqV/MXTMHQxmKp/fcrZjRUVQtMREQ08FQVFuM/hXqpnWzMx5SZEwUmInJdLKiIqF8plUr877QYKKwWAECxxhcfrz8kOBUR0cBhs9nw7paTqFHZCyqtxYgfzU+AXM6PhUTdwX85RNTvBg+NwX0ejgkpNspjcPkYn6ciIuoPJ3cfxn7NIKn9aFA9gsOCxQUicnEsqIhIiIcXpiLMVAUAsMrk+HNWLYy1dYJTERG5t7qSUqzMVUrtIaYyLJzLoX5EPcGCioiE0GpU+HFykNTO1QXhsy/3C0xEROTebDYb3t+YgVKNDwBAaTXjJ2mxUCr4cZCoJ/gviIiEGZ0Qi/nqEqn9hS0G1zPPCkxEROS+zqQfw3ZNrNR+yLcGgweFdXAGEXUGCyoiEuqxRSkIMtln+bPIFXjnVAVM9fWCUxERuZeGigqsvGaT2jGmCjx4d4rARETugwUVEQml16nxv+P8pPYNXQjWcegfEVGv+mD9ERRp7L9r5TYLnp4WBbWSHwOJegP/JRGRcOPGDMVspWPo36eWKOSeOS8wERGR+7hwMAOblYOk9v2eVYgbEikuEJGbYUFFRE7hu4smwN9UCwAwy5V4J6MYZoNBcCoiItdmqK7Gn7NNsMnsH/kiTZX45kLO6kfUm1hQEZFT8PTQ4seJXlL7ij4MG7/cJy4QEZEb+HjdQeRrAwAAcpsVT08Kh0alvMNZRNQVLKiIyGlMGD8MaXLH0L+PjBHIO3dJYCIiIteVfTQTGxSDpPY9+goMGxYtLhCRm2JBRURO5fv3jIev2b7Ar1Ghwp+P5MPcaBSciojItTTW1OLP5xpg/WqoX5ipCo8u4lA/or7AgoqInIq3px7/M0IntS/qw7FhPWf9IyLqig+/PIhbukCp/dSEYGjVKoGJiNwXCyoicjqTJo7ENFmx1P6wMQw3L1wVmIiIyHWcyziDDTLH0L6F6lIkjBosMBGRe2NBRURO6Yf3jIdfs1n/3jpSACOH/hERdai+phbvnG+QZvWLMFZg6T0c6kfUl1hQEZFT8vbywFMjNFL7hjYIazYcFpiIiMj5/XdDBorUvgDsC/g+MzEYWq1abCgiN8eCioicVnJKIuba8qT2F41BuHzuisBERETO68TRLOxAuNR+UF2MYaOGCExENDCwoCIip/bde1MQ0lgBALDKFHg7owSG+gbBqYiInEt1eSX+csmxGPpgQwkevneKwEREAwcLKiJyanovTzyd6AmZzQoAyNP4Y/Xag4JTERE5D5vNhn+sP44KlX1xdKXVjGWTw6DWcKgfUX9gQUVETi9hQiIWKwul9iZZFM4eOiEwERGR8zi4/RAOqqOk9qNeZRg0Ik5gIqKBhQUVEbmEb98/BVHGcqn9f9lm1JWUCExERCRe2c1c/L1QL7VHGItx76KpAhMRDTwsqIjIJWg0GiybEgGF1QIAKNX44l/rj8NmtQpORkQkhtVkxF92XEKtyl5QaS1GPD13OJRKheBkRAMLCyoichlx8TF42LdGau/VxeLApr0CExERibNp7R6c1DsW8P1upBnhEcECExENTCyoiMilPLRgIoaZy6T23yr8UXz5msBERET978apLKwyRkjtZFsJ5s0cJzAR0cDFgoqIXIpSIcdP54+AztIIAKhX6vDW/hyYDYY7nElE5B4M1dVYkVkJk1wFAPA11+Gpe5Igk8kEJyMamFhQEZHLCQvxxw9jHR8cLugj8MUXewQmIiLqHzabDe9/fgC5uiBp39NjvOHn4yEwFdHAxoKKiFzSzCmJmCYvldqf2GKQfSRTYCIior53YtdBbNEMkdqLtGUYP3aowERExIKKiFySTCbD/yxORpCpGgBglSmw4kIj6svK7nAmEZFrqsi9jT/f1kntGFMFli5OEZiIiAAWVETkwjw9tPhpSjDkNvvU6YVaf/zzy2OcSp2I3I7VZMI72y+gSu0JAFBbTXjurlhoVErByYiIBRURubRRIwbhQe8qqb1XG4uDW/aJC0RE1Ac2r2s5RfpjYUbERIcITERETZzqa43s7GysXbsWly9fhtlsRmRkJObNm4e0tLQuXae6uhrHjx/H1atXcfXqVdy6dQtWqxXPPPMMpkyZ0uG5t2/fxpo1a3D+/HkYDAaEhoZi5syZWLBgAeRy1p9EzuiRhSk488ERXFYGAABWlvkh/up1hMTFCk5GRNRzN06dw/uN4dLX4ONtJVg4a6rYUEQkcZoKISMjA8uXL8fp06cRExODsWPHorCwECtXrsT777/fpWtdunQJ//jHP7B7927cvHkT1k4O/7l8+TJ+9atf4ejRowgJCUFycjJqamqwatUqvPXWW7DZbN15a0TUx5QKOZ6dN7zFVOp/2ncTJk6lTkQurqGqGm+crGoxRfpPFo3jFOlETsQpeqhqa2uxcuVKWK1WPPfcc0hJsT9gWVlZiZdeegmbN2/G+PHjkZCQ0Knr+fr6Yu7cuYiLi8OQIUOwfv16pKend3iOxWLBn//8ZzQ2NmLp0qVYtGgRAMBgMOB3v/sdjh07hn379mHmzJk9e7NE1CfCQgPwo8F5eDvX3r6sC8OHX6Tj8W/PFRuMiKibbDYb/rHuKPK0jqF+T4/2gp+vp8BURPR1TtFDtWfPHtTX1yM5OVkqpgB7YbRkyRIAwKZNmzp9vfj4eDzxxBNIS0tDVFRUp77FycjIQFFREWJiYqRiCgC0Wi2+//3vdzkDEfW/mdNG4y5ZkdReh2gcP3hSYCIiou7bveMo9qocxdT92hKMHxcvMBERtcUpCqrMTPvaMampqa1eS0pKgkqlQlZWFoxGo5AMgwcPRkhICG7duoXi4uI+y0BEPffD+1MR1Vgutf/vGlBSWNrBGUREzif3+i38s8ixWO8wYzG+fe8kgYmIqD1OUVDl5trH6MTGtn6AXKlUIjo6GiaTCfn5+X2W4ebNmwDsxVNbmvY3HUdEzkmn0+Dn08Khtti/gKlR6vHm1vMwm8yCkxERdY6hrgFv7MtFo0INAPA0N+Bnc+OhUjrFkxpE9DXC/2XW19ejrq4OAODv79/mMf7+/rh27RpKS0sxaNCgPslRWmr/BjsgIKDdDM2Pu5Nnn3221T61Wo3XXnsNABAYGNidmL1G+dUv5aCgIKE5yDm5+v0RFBSEn+RuwZuF9g8jF9Uh+Gz9ITz95EOCk7k+V783qO/w3ugdNpsNr3y4BrmacGnfL0eoMGpsosBUPcN7gzriDveH8B4qQ7NZuDQaTZvHNO039OGMXU3Xbi+DVqvt8wxE1Hvuf3g+ZlodvdqfGoJwYOs+cYGIiDph02dbsR2OYuoBVSHuWsgJsYicWa/0UL355pu4detWl8556qmnEBcX1xs/vl90dcr0FStWdPh6aWmp0GnYm74FKCkpEZaBnJe73B8/vD8V2Z+dQb7aDzaZHL8/V4u3gzIREBN955OpTe5yb1Dv473Rc7cvXMZbuWrp09lQYwm+9dAUl/8z5b1BHRF9f8hkMoSFhfXoGr1SUJWUlHT5+abGRvt6MU09P0379Hp9p47tbVqtFnV1ddLPEpGBiHqX3lOPn0+NwC+OVsMkV6Ja5YkVOy7j5SVBUOp0ouMREUkMlVV443ABDLoQAICHuQE/mzccao1KcDIiupNeKaiangvqDr1eD71ej/r6epSXl7dZUJWX22fs6svnjgIDA1FXV4eysjLExMQIyUBEvS92aDSeyM3E3wrtv+7O6SPx4ae7sfSxhVwYk4icgtVixj8+P4gc3RBp30+GqxEa7rrPlBANJMKfoQIgFTDXr19v9ZrZbEZubi5UKhXCw8Nbvd7bGW7cuNHm60372yq2iMi5zbsrCVPljgll1qricGTLPnGBiIia2bF2F/Y0K6YW6soxKWWUwERE1BVOUVAlJSUBAI4ePdrqtZMnT8JkMiEhIQFqtVpIhhs3bqCoqAiRkZEIDg7uswxE1DdkMhn+9/4URJoqpX3vlPrhdtYFcaGIiABcPnIC/zJESu1hljI8vrj1mphE5LycoqCaNWsWdDodTpw4gWPHjkn7q6qq8MEHHwAAFi1a1Oq8ZcuWYdmyZdJwvJ6YOHEigoODcfPmTWzatEnabzAY8O6777abgYhcg16rwi/nDIHOYn8eskGpxWvHytDQC78/iIi6o+p2Hl6/aIFZbh+S7Guuw/P3JEKtdIqPZ0TUScLXoQIAT09PPPnkk3jrrbewYsUKjBw5El5eXsjKykJdXR3uvvtuJCa2Xn+haSIMs7n1gp2//vWvpe3CwkIAwKeffootW7YAsC/U+8QTT0jHKJVK/OQnP8Fvf/tbrFq1CkeOHEFgYCAuXbqEiooKTJgwAWlpab35tomon0VFBOHp4aV4/Yq9fUsXhL+uPYZnH58DORfMJKJ+ZDYY8ObWCyjVRwEA5DYLfpbsjwA/b8HJiKirnOYTRGpqKl5++WWsXbsWV65cgdlsRkREBObNm4eZM7u+/sKVK1da7SssLJSKK5Wq9aw5w4YNw6uvvoo1a9bgwoULyMnJQUhICBYtWoSFCxdCLuc3RkSubvLEEbg//wjW1fkBAA7oBiP+851Y/MjdgpMR0UBhs9nw8ae7cUbveG5qaVADEkfxuSkiVySziVwMaQArKCjgOlTktNz9/jBbrFj+4WGcU9hn7VRYLfjt4DqMmjpRcDLn5+73BnUf743OO7ptP14tC5Hak1CC5x+d6rYzj/LeoI6Ivj96Yx0qdrkQ0YCjVMjxs3vHwd9UCwCwyBV44wpQfjNXcDIicnf5Fy7j/4p8pHaEqRJPP5DitsUU0UDAgoqIBiQ/Hw88PykYSqsFAFCh9sYfd1yFsaFBcDIiclf1FZV47XAh6pVaAIDW0ohfzoqFXtd3sxgTUd9jQUVEA9bwYdH4frhBal/Uh+Nfn6bDarUKTEVE7shiMePttcdxU+dYfuXpeCWio7gcC5GrY0FFRAPa3Xcl4S55sdTeoYrB1s0HBSYiInf06RfpOKaNktr36cowJZWTUBC5AxZURDSgyWQy/M+DkxBvdDwM+25lAM5mctFfIuodh9JP4lNTuNROMhfiO/dOEpiIiHoTCyoiGvA0ahV+uXAE/I3VAOyTVPzxnAGFN/MEJyMiV3f90jX8X45jlZqIxnI890AylAp+BCNyF/zXTEQEICA4EL8a5wWV1QQAqFHq8YfdN1BfXS04GRG5qsqiUvzhSCkaFfZJJ/TmBrwwPRyeXp6CkxFRb2JBRUT0lfixI/C/IbVS+6YmEG9/dhQWk0lgKiJyRcaGBry+8SxK1PYp0uU2K34eZ0Nk3CCxwYio17GgIiJqZubcSbhfVSC1j2mj8cknO4UuxE1ErsVmteDfH+/FBZ3juaml3uVImpIkMBUR9RUWVEREX7PkwRkYby6S2mvksTi4YZfARETkSrau2YbtmlipnYYi3HfPFIGJiKgvsaAiIvoapUKOZx9KQYSpUtr3TlUoLh08Ji4UEbmEE9v341/mQVJ7qKkMP354MmQymbhQRNSnWFAREbXB00OLF+YNhafZvvCvUaHC768oUHgxW3AyInJW14+fwhuFvrDKFAAAf3MtfnVvIjQqleBkRNSXWFAREbUjMiwAz0/whdJqAQBUqz3xyuES1BQWCk5GRM6m7HoOfpdlhEGpAQBoLUb8ZkYkAvy8BScjor7GgoqIqAOjRw7C/w62Su08bSBe35gFY21tB2cR0UBSX1aG3+65iTKNY0a/nyXqMGRQqOBkRNQfWFAREd3BXVMT8bCPYz2qLH0UVn6SDiunUyca8MwGA95cewI3dCHSvifCGjFh3FCBqYioP7GgIiLqhEcXTsB0RanU3quLxWef7OB06kQDmM1qxX8+2o0T+hhp3z36ciycNU5gKiLqbyyoiIg6QSaT4ScPpmKk2VFUfSQfgvSNewSmIiKRNq7Zgc2aIVJ7Ikrx+OJUgYmISAQWVEREnaRWKfHLB5IQZqqS9r1TFYJzB48LTEVEIhzZno7/mqOl9hBTOZ59KAVKBT9aEQ00/FdPRNQFPl56vDRvCLzM9QAAs1yJ319T4sb5y4KTEVF/OZdxBiuK/GCV2T9GBZpq8Jt7E6HTcHp0ooGIBRURUReFhwXihQl+UFvtk1LUK3V4+XgVCm9xOnUid3cjOwd/uGiFUWEvnvRmA15Mi4S/n5fgZEQkCgsqIqJuGDlyMJ4bYoPcZp9SvULlhZd33UBlWYXgZETUV4ryCvHKkRLUKXUAAJXVhF+P1WPQoDDByYhIJBZURETdlDp5NP4nwFFA5av98Nv1Waiv4RpVRO6mqqQML2+/jnKVvSdKbrPiuWgjEsbEC05GRKKxoCIi6oF5d0/Bo6o8qX1VE4zX1xyDscEgMBUR9aaGqir8dsM55Gn8pX0/8inBpBnjBaYiImfBgoqIqIe+8dBMLJDlS+3T2gj8+aN9sBiNAlMRUW8w1dXh9c8ycEXrWLj3UW0h5t8zQ2AqInImLKiIiHpILpfj+w/PwBSrY1KKdO0g/Hf1dlhNJoHJiKgnLAYD/vzRPpzSRUn77lYU4hv3TxeYioicDQsqIqJeoFQqsOyRKRhtKZH2bdQOxacfbIXNahGYjIi6w2o04p8f7MR+vWPh3skoxhPfmA65nB+fiMiBvxGIiHqJWqXCL7+RglizY6KKT9TxWPvhFthsNoHJiKgrrGYz3vtgO7bphkr7Eqxl+OnDk7lwLxG1wt8KRES9yEOnxksPjEWkuUrat0o+FJs/YVFF5ApsVis++XAr1mscxVS8pRy//sZEqFVKgcmIyFmxoCIi6mV+Xjq8sngkQk3V0r5/WYdg19rtAlMR0Z3YbDas/WgLPlU6iqnB5gq89FAS9FqVwGRE5MxYUBER9YEAPy+8sigegeYaad9fG6Kxf8NugamIqD02mw2bP92GVbI4aV+kqRL/74Ex8NJrBSYjImfHgoqIqI+EBPrilbmx8DXXAQBsMjnerg7Fka37xAYjolZ2rt2Jf1kGS+1QUzVeuXcUfL30AlMRkStgQUVE1IciwgLwyl2R8DI3AACsMgX+VBqIE7sPC05GRE32b9yDlQ2RUjvQVINXFsYjwM9LYCoichUsqIiI+lhMVAhenhYMvcUAADDLlXgt3wuZ6ccFJyOi/VsP4O2qENhk9o9EfqZa/HbeYIQE+YoNRkQugwUVEVE/GBIbgeUpftBaGgEAJrkKf7ipw/HDp8UGIxrA9uw8irfLAmCVKQAA3uZ6vHJXFMLDAgUnIyJXwoKKiKifDB8Wg5eSPKWiyixX4rVrShw7eFJwMqKBZ9eOo3inyBvWr3qmvM31eGVaCKKjQwQnIyJXw4KKiKgfjUoYgv83RgOd2TH87/UbGhzec1RwMqKBY/umdPyl2Fsa5udjqsVvJwdgcGyE4GRE5IpYUBER9bMRY4bj5UQV9F9NVGGRK/BGvhcObE0XnIzIvdlsNmxduwsrq4KlYsrXVIvfTfLDoKExgtMRkatiQUVEJMCwpFF4ZbwnPJrN/reiLAD7v9wpOBmRe7KvM7UVf282m5+/qQa/nxaM6GFDBCYjIlenFB2AusZms/XqdXrreuRe+ur+kMlkvXo9Vzc0YSheUeVg+dEK1Cp1sMoUeLs2AuY1m3HXNxbwz4uol9isFqz/YBP+qxgm7Qsw1eB3d0UhPDpUYDIicgcsqFyA2WxGQ0MDjEZjr33AraioAABYrdZeuR65l766P2QyGdRqNXQ6HZRK/voBgLhhg/BblQovHShGjVIHq0yOd0xDUPvhRix+dCFkcoXoiEQuzWoy4sPVW/G5xlFMBZpq8Lt5gxHG2fyIqBfwE42TM5vNqKqqglarha+vL+Ty3hml2fRh1mw298r1yL301f1htVphMBhQVVUFHx8fFlVfiY2NwO9Uary05xaqlHoAwH9k8ah5fzMeXXI35CqV4IRErslsaMA/P9iF7TpHMRVsrsFvFw5FKNeZIqJewk8zTq6hoQFarRYeHh69et2moUQcUkRt6av7Q6FQSPdyQ0MDvLy8evX6rmxQVBBena/E8m1XUaK0/7l8po5H9art+OG3Z0Op1QpOSORajLU1ePvjAzikHyrtizJX4uV7ExHg27v/TyWigY2TUjg5o9EILT9IkZvRaDQwGo2iYzidiBA/vHbvCESZK6V927VxWPHBXhjr6sQFI3IxDeXl+P3HR3BIHyvti7eU4w8PjWUxRUS9jgWVE7PZbLDZbL02zI/IWSgUCun+ppYCfT3x+wfHYqi5XNp3SDcYv//oMBq+eraNiNpXXVCEl744g9P6aGnfWFsZXnlkArw9+AUlEfU+flInInIyPp5avPLNCRhjLZX2ndZHYfnnp1FVWCwwGZFzK83Jxa83X8ZlfZi0b4q8FL9+ZBJ0aj6LSER9gwUVEZET0mtV+M0jqZgscxRV2fow/HLLVeTfuC0wGZFzun7hKn6+txC5uiBp31x1GZ59eDLUSn7cIaK+w98wREROSq1S4rlvTsZctaOoytf44/n9xbh49rLAZETOJfPoWfzqRB3K1d7Svoe8KvDjhyZDqeBHHSLqW/wtQ0TkxJQKOX780BR8y7NM2let0uPFM404sPuowGRE4tlsNmzblI7fXVXAoNAAAOQ2K34YWIXvLJ7EmWyJqF+woCKXVV9fj3/+85946KGHMGbMGAwaNAgjR47EPffcgzfeeAN5eXmiIw4YKSkpiIiIEB3DbclkMjxy7xQ8E1QBpdW+NphJrsKfCn3xxSfbYbVwPTkaeCwmI97/YDv+VhUMq8y+ALbW0ogXBhmwcF6K4HRENJCwoCKXlJmZialTp+Lll1/G6dOnMWzYMCxcuBDjx4/HzZs38fbbb2PatGlIT08XHbXfHT58GBEREVi2bJnoKNTL7po7Cf9vuA2e5gZp3ypLDP723naYa2sFJiPqX4aKCvzpvd1YJx8k7fMz1eAPSTpMmJokLhgRDUhc2JdczoULF/Dwww/DYDDgf//3f7Fs2TLo9XrpdavVim3btuH3v/89CgoKBCYl6n2JExLxWlAefrv3NorUPgCAHdohKP7wEH62YBS8oiIFJyTqWxVXr+LVPTeR7TFY2hdjLMdv7o5HcGigwGRENFCxh4pcis1mw9NPPw2DwYDnnnsOL7zwQotiCgDkcjkWLFiArVu3YsyYMYKSEvWdqEER+OMDIxFvdjxXddozBr/YnoPczNPighH1sewDR/FsehmyPRxDjMdZS/HqN8ezmCIiYVhQkUvZt28fLl68iLCwMDz99NMdHuvt7Y3hw4dL7YaGBrz11lu46667MGTIEAwfPhwPPPAA1q9f3+b5zZ8Leu+996TzUlNTsXLlSmlR2qysLCxduhSjRo1CfHw8vve97+H27dbTWi9btgwRERE4fPgw9uzZg/vuuw9Dhw7FyJEj8cQTT+Dq1autznnzzTcRERGBTz/99I4Zm37GN77xDQDAZ599hoiICOm/N998s8W5t27dwi9+8QukpKRg8ODBSExMxA9+8ANcuHChzZ9lNpvx5z//GVOmTEFsbCwmTZqEP/7xjzAajW0eT33L18cLv/1WCiYpHIv95usC8YtzwNFNu7hoMrkVm9WKXZ9txa9zPFCu8ZH2z9VV4tePToaHXiMwHRENdBzyRy5l9+7dAIBFixZBqez87VtbW4tvfOMbOHv2LAICAjBr1iw0NDTg0KFDOHbsGDIzM/HKK6+0ee7y5cvxwQcfYNy4cYiKisLRo0fx+9//HvX19ZgxYwa+9a1vISoqCpMnT8aFCxewfft2ZGdnY9euXdDpdK2ut2nTJqxatQpjxozBnDlzcPHiRWzduhWHDh3C559/jlGjRnXvDwfAxIkTUVJSgn379mHQoEGYMGGC9Frz62ZkZGDp0qWoqanBsGHDMGfOHBQWFmLr1q3Ys2cPPvzwQ0ydOrXFtX/84x9j8+bN8PDwQFpaGmw2G/75z3/i3Llz/PAuiFatxC++mYo120/i4zIPAECDUotXqyLxyKrNePhbc6FQqwWnJOoZU30d/vvJXmzWxElfA8ttVnw/yoyF01M4kx8RCceCykXZbDagoa77539VjNjM/Tw7mM6jR//zO3fuHAAgMTGxS+e99tprOHv2LKZNm4Z3330XHh72D59Xr17Fgw8+iHfffRczZszArFmzWp27adMmbNmyBcOGDZPOmTt3Lv7+97/j888/x/PPP48f/OAHAACj0YglS5bg0KFD2LBhA775zW+2ut7777+PP/7xj/j2t78NwP53+eqrr+Kvf/0rnnvuOWzbtq1L7625Rx99FIMGDcK+ffswYcIEvP32262OqampwY9+9CMYDAb84x//wKJFi6TX0tPT8dhjj+Gpp55CRkYG5HL7p5cvv/wSmzdvRkxMDL744guEhYUBAHJzc/HAAw/wWTWB5DIZHpk/HoNPXcVbWXVo+Grq6E+Ucbj+/l4su288PII4FIpcU2VePt7YcgHn9HHSPm9zA34+MQCjR0QLTEZE5MCCylU11MH6zKPdPl3UIC35/30E6D27fX5FhX14U0BAQKfPqa+vx8cffwy5XI4//OEPUjEFAHFxcXjmmWfw4osv4j//+U+bBdUvfvELqZhqOmfWrFnYsmULIiIipGIKANRqNZ544gkcOnQIR44cabOgSk5OloopwD4l9s9//nOsW7cOWVlZOHHiBJKTkzv9/rrqk08+QXFxMZ566qkWxRQATJ8+HUuXLsW///1v7Ny5E/PmzQMArFq1CgDw85//XCqmACA6OhrLli3D888/32d5qXNSxsXhj8Gl+MOuayhQ2odEZehj8IsNl/CrScGIHBkvOCFR11zLPItXz9SjRO+YaGWwuQK/WjQKIQHeHZxJRNS/+AwVuZTuDC07e/YsDAYDxo4di9jY2FavP/jggwCA48ePt3n9adOmtdoXHR3d7msxMTEAgOLi4jbz3Hvvva32qVQqLFiwQMrRl5qmkp8/f36br0+cOBEAcOrUKQCAyWTCqVOnIJfLsXDhwlbH33fffX0TlLosOiIQb3xjDMbbSqV9t7WB+NnxehzYxUWAyTVYrVZs3ZiO5y/IUKLxlfZPV5TitW8ls5giIqfDHipyKf7+/rh27RrKysrufPBXioqKAABRUVFtvu7j4wNvb29UV1ejpqYG3t4t/2fdvEemSdPMgh291t5kDZGRbU9r3ZSvsLCwzdd7S9OEGV/vnfq68vJyAPZeQaPRiJCQEKjbeB7H09MTPj4+qKqq6v2w1GVeei1e+NZkfLzpGD6v9QNgf67qT0VanP1wD77/wCRo23i2j8gZ1FXXYuWXGTioCG/xvNTSwDrcN28Kn5ciIqfEgspV6Tzsw+e6qWlCB7OAZ6h6YtSoUTh+/DiysrKknqXO6sz/iNs6pr/+B96d3jer1drlcywWCwB7QdXWpBmAfer5pKSkFrn4QcZ1KBVyfOfeSYg9nIW/XLGgXqkFAOxAOLI/ysQvpoQicnjcHa5C1L+unjqPP52qRoEmXNrnbarDs2O9MG7shA7OJCISiwWVi5LJZD16Fkn2VUEl6++CqodmzZqF9957D5s2bcJvfvObTs30FxISAsA+gUJbqqurUV1dDb1eD0/P7v+ZdlZbU6oDQF5eHgAgNDRU2qdSqQDYnwP7OovFgpKSki7//LCwMFy7dg3PPPMMRo4c2eYxzQtuf39/qNVqFBcXw2g0tuqlqq2tZe+Uk5oyORGDI/Lwp13XcE0bDAC4qQ3Ecxl1+J9Lu5C2eCZkcoXglDTQWU1GbFm3B/9tjIRZ4yftTzAU4NlFiQgI4aQqROTc+AwVuZSZM2di2LBhKCgowDvvvNPhsTU1NcjOzsbo0aOh1Wpx+vRpXL9+vdVxa9euBWB/dqg/emE2bNjQap/ZbMaWLVsAoMWEFE3FYFu5Dx06BJPJ1Gp/UxHW1BP1dU3PfW3fvr1TeVUqFcaOHQur1SplbK69dbzIOYTHROC1b6dgkTxf2mdQaPB2XST+/O/NMBQXCUxHA11Nbi5e/88u/Ms0CGb5V1/02az4pqYQLy+dxmKKiFwCCypyKTKZDO+88w60Wi3efPNNvPrqq616b2w2G3bs2IG7774bp0+fhl6vxyOPPAKr1Ypf//rXLY6/du0a/u///g8A8N3vfrdf3sPx48fxySeftMj7pz/9CXl5eRg5cmSLtaNSU1MB2Iu+W7duSftv3ryJ3/zmN21ev6mH69q1a22+vmTJEgQEBODPf/4zPv3001ZDDevr67FmzRrk5+e3OAcA3njjDemZNMDe29bW1OzkXNRaDX7wrbvwy1gTPMwGaf9uj3j8dMMVXNp7kGuJUb+y2Ww4tXU3ntlVgKOejsmCfE11eDlBgUcfSoNSxUE0ROQa+NuKXE5CQgI++eQT/OAHP8Bf/vIXvPvuuxg/fjyCgoJQXV2Ns2fPoqSkBFqtFuHh9rH4v/rVr3Dy5Emkp6dj0qRJSE1NRX19PQ4fPgyDwYDvf//7mD17dr/kX7p0KX72s5/hgw8+QExMDC5evIjs7Gx4enrirbfeanFsTEwMHnroIXz++eeYO3cuUlJSUF9fj5MnT2LWrFlobGxsNYQwKioKI0aMwJkzZ7Bw4ULEx8dDoVBg7ty5mDt3Lnx9ffHuu+/i8ccfx7PPPosVK1Zg2LBh0Gg0yMvLw5UrV1BfX4/du3cjONg+TOyBBx7A1q1bsXXrVkyfPh1Tp06FzWbDgQMHkJqaCplMJg1ZJOc1aVIiYuMq8aetF3FZYR9ala8LxK/yrHjwvQ14+KGZUHtxBjXqW4ayMry37hC26uIBjWP/GFs5fnr/GPj59OxZWyKi/sYeKnJJEyZMwKFDh/Diiy9i7NixuHjxIjZu3IgTJ04gMjISzz77LA4cOCANb/P09MQXX3yBn/3sZ/D398fOnTuRkZGB0aNH469//SteeeWVfst+zz334L///S/kcjm2b9+OgoICzJs3Dxs3bkRCQkKr49944w089dRT8PT0xP79+5GXl4ef/OQnWLlyZbs/41//+hfmz5+Pmzdv4vPPP8fHH3+MrKws6fUJEyZg9+7d+NGPfgStVotDhw5h//79qKmpwezZs/HPf/4T8fGOdYtkMhn+9re/4fnnn0dAQAD27NmD8+fP43vf+x7+/e9/c8IKFxIS5Is/fDsF9/vWQvZVr5RVJsdn6mH4xScnkZN5RnBCcmfZh47hp19eshdTX1FaLVgS0oj/9+1JLKaIyCXJbE40ziM7Oxtr167F5cuXYTabERkZiXnz5iEtLa1L16mursbx48dx9epVXL16Fbdu3YLVasUzzzyDKVOmtHlOSUkJTpw4gdOnTyMvLw/l5eXQ6XSIjY3FvHnzen2h1YKCgjsOsbHZbCgrK0NAQECvf2AVNsvfALZs2TJ89tln+OyzzzB58mTRcTrU1/dHX97b1Hnnr+bh/w7lo0jpJe1TWs14VHkb9z4wE0qNptU5QUFBANCtCVHIvXV0bxjr6/Dp5/uwVj4YVpnju9xoSxWWzYzFkKigfstJ/Y+/N6gjou8PmUzW5jI4XeE0Q/4yMjKwYsUK2Gw2jBgxAl5eXjh37hxWrlyJmzdv4rHHHuv0tS5duoR//OMfXfr577zzDrKzs6FWqzF06FDExcWhqKgIZ86ckYZOdSUDEZGzGxUXgbejgvHexuPY3ugPADDLlVhlHYTjqw/hqWkxiBw2RHBKcnXXzlzEn4+X4IbOcS/JbFbc512DRxckQ63kTJNE5NqcoqCqra3FypUrYbVa8dxzzyElJQUAUFlZiZdeegmbN2/G+PHj2xwO1RZfX1/MnTsXcXFxGDJkCNavX4/09PQOzwkMDMT06dMxbdo0aLVaaf/JkyfxxhtvYPPmzRg7dizGjBnT/TdKRORk9BoVfvzQZKScvoq/nK5Eucq+dMBFXTieyWjAQ+f24cGFk6HWtl7Umagj9bV1+HhTBjaZQ2DVBUv7Q01VeGZSGEYOa3vZBiIiV+MUz1Dt2bMH9fX1SE5OloopwF4YNc0utmnTpk5fLz4+Hk888QTS0tIQFRXVqSFFzzzzDObMmdOimAKApKQkzJw5E4B9mmoiInc0fmwc/u+hkZgmK5b2meVKfGIIxbKPjuPc4UzOBEidYrPZcGzPEfzks3PYYAlrMcRvnqoEb31zHEYOixaYkIiodzlFQZWZmQnAMUV0c0lJSVCpVMjKyoLRaOzvaADsM60BQEVFhZCfT+7h7bffRl5entM/P0UDl7enHj97dDp+ObgR/sZqaX+eJgC/vuGBP/9nK8raWBONqEnBhQt4/d/b8YcCP5SqfaT9YY0VeHkE8OOHp0H//9u78/ioqvv/46+ZzGSZJCQsSSAgBAIaWWSPCAouICKoaGXR/myxpT7QKiKK5VGs1YdaLBSwaqv9qigFN6xoq1JUEJCALAookhCzQjCBLGQSSEgy2+8PmkhMgmFyJ5OB9/Mf4d47Z87Vj3fO5557PydMs50icm5pEwnVoUOHAOjVq1eDfRaLhe7du+NwOOqti9OaatfdiYqK+okjRUQC32UjB/L8zUlM5HtMHnfd9g2hvfj5u+ms+b9VuCor/NhDaWucZXZWPvca/2/t93xhS6jbbnE7mWLJ55lplzBoSJL/Oigi4kN+f4eqsrKSiopTP8wdOnRo9JgOHTqQlZVFcXExCQkJrdg7qKioqHv/6vQFV3/K3LlzG2wLDg7m6aefBk69s/VTPB4PpaWlWCwWwyuh1bZXW81N5HS+jg+Px4PZbCYmJkZV/tqomJgYFsxJ5Ia9aSxen0GW5dS6VeXWCP5SGUGfVTu4d0AUwyZdi8ncJu7NiR94XE4+f/cj/p7pJM+WWG9U0c9Vwu8mXkLvpKv810FpE2p/S2qruYmc7lyID7+Ppquqqur+HNJIid7Tt59+bGt56aWXKC8vp0+fPiQnJ7f694uI+NMlgy5m+YCLeOv9z1l+0E212QpARng892fD5U+v4rfXD6LHoEv83FNpbd9t28FzmzL5ytYdbD9sD3dVcffF4UyecBNm3TARkfOAIQnVkiVLyMvLO6vP3HvvvfTu3duIr/eZ999/n23bthEREcHs2bPP6k760qVLz7i/uLi4WetQud1unE6n1qGSVtUa61C53W6Kioo0QxUgrhvdj8GlJ/jnZ2mkVP2wblWKrSc7Pivmhg0ruPX6Swnv1NGPvZTWUJZfwJsf7+ZjawJu2w/FJUweD+MiK7n9mn60jwilpLjYj72UtsTf6wxJ2+bv+Ggz61AVFRWd9ftN1dXVAPWq6lVXV2Oz2Zp1rK9t2rSJN998k5CQEObPn09cXFyrfbeISFsU1z6CRbOu56tvc3jm42/JCooGwGG2soZebPggh+ntv2HshFEEh6jwwLmmquIka9du5V+VHakIrr8+WV9PKXNvHErfxK4aNIvIeceQhKr2vSBv2Gw2bDYblZWVHDt2rNGE6tixY0Dz3jsywq5du3jxxRcJCgrioYce4sILL2yV7xURCQRD+/fkL7HhbNy2n1WZVRyznFq7qiw4gn9URLDmjT1MiXNw9TXDsTbxKLcEjqqKCj5e/yVrSsOxW7vVGznEOcqZ0T+ay4aOIDY2tulGRETOYX5/hwpOlSVPS0sjOzubbt261dvndDo5dOgQVquV+Ph4n/dl//79PPPMMwDMnj1bC/mKiDTCbDJxzaj+jBzmYM26L3m/LIKaoFPvVxUFR/H3Unj39d1M6VDJVWMvxRIR4ecey9mqttv59NOd/OtEe0qD48D6wz6bq4opMTVMGjuEYGubGEqIiPhNmyjNNGTIEAC2b9/eYN/u3btxOBz079+f4GDfPkKSnZ3NokWLcDqdzJo1q9F1sURE5AdhIVZ+ftNl/H1CN67le4Lcrrp9R0Pa83xFV+59+2s2rP4IZ+kxP/ZUmqum6Chr3/iAu9ek81JNd0qDf3hnLtjl4IagAl64KZFbJiQrmRIRoY0kVNdccw1hYWF8+eWX7Nixo257WVkZq1atAmDSpEkNPjdnzhzmzJlT90hgS+Tn5/OnP/2JkydPMmPGDK688soWtykicr6IiWnPb39+DX8fHc1YCjB7fkisCkI78qwjkbvXHODfK/9DRf73fuypNKUsJ4fVr/2buz44yD88fSgJ+WHtRavbwURzAS+Oi2Pm9KuIjoo8Q0siIucXk+enSs21ku3bt7Ns2TIA+vbtS2RkJPv27aOiooIJEyZw5513NvjM1KlTAXj++ecbPLu9YMGCuj8fOXKE48eP07lzZyIjT/0I9OzZk5kzZ9Yd8/DDD5Obm0u7du0YPHhwo33s2rUrkydPbtF51iooKGhWlb+SkhI6duyoKn/Sqlqjyp+vYlt8q7nVmPILili9OZ3Nzg64TfXv3dmcVYxzH2bSyCRiL2rb1V7PB9/vS+U/O7P5zNqdmqD6T4JY3E7Ghdq59aoBdOp05sXt/V2pS9ouxYacib/jo81U+TPCiBEjePzxx1mzZg0ZGRk4nU66du3K+PHjueqqs18UMCMjo8G2I0eOcOTIEQCsVmu9fbWLC5eXl7N58+ZG2+zbt69hCZV4r2vXrj95zJQpU+rehROR1hffJYY502O4tbCM1Z8fYEtVZF1iVWkJ5d/05oNd1YxKWceNA+PpM3SAkutW5Ha5SN3xNf9OLWJX6AV4Qusntla3k6siKplyVT9io8P91EsRkcDQZmaozjeaofJebUI1ZcqUJo9JTk7m9ttvb60unbUlS5awdOlSli5dyrRp0/zdnQY0QyVN8fZOYmHZST7cvI9P7SFUBjWs/NezuphxMW5Gj+xHZEetZeUrpQVH2bT9AOtLrRwO6dBgfzvXSSZ0cjFhzADah59dhUZ/32WWtkuxIWfi7/g4p2aoRM6WZqBEAkdsVBi/ujGZ6dUO1n/+DR9876LQ2q5uf05IJ/6vHF77KJ8Rjl2MSwin/2VDMNs0O9JSzuNl7Nm6l/WHT7Ir5AJc5jj4Ua7U1WHnpoQwxoy6hFBrkH86KiISoJRQiYhIq7GFWLlx3FCud7nZsWM//063kx4cU7e/JsjK50EJfH4EOr+xlytNhVzRtwtdhw3BZNViwc3lPllJ7q49bP3uKBvN8ZSExEFYw+MucRRyY/9Yhg69FLNmi0VEvKJH/vxEj/x5r/aRv++/b16lsCeffJIXXniBG264gRdffLHevuLiYsaOHcuxY8dYs2YNw4YNA+Do0aO8++67bNiwgdzcXEpKSoiOjmbYsGHce++9DBo0qNHvqqys5JVXXuHDDz8kJycHj8dD165dGT16NHfddRfdunXj0ksv5fDhw41+/p133mHkyJHN/DfhO3rkT5rii0czDqZnsX7PQTZVR1Nuabi4O0D3yqOMDDnOyH7duGBgf8wW3Q/8Mbejhqyv9rHtQAFfONtTENb4o5MdHCe4ynacscN6E9/rAsO+39+P7UjbpdiQM/F3fOiRv/OYx+OhwuH2+vMW96lBrNPp+okjjRVuNbf6APrhhx9my5YtfPDBB1xzzTX13r168MEHKSoqYu7cuXXJFMDHH3/MU089RUJCAklJSURERJCbm8t///tf1q9fz4oVKxgzZky97zl69CjTp0/nu+++Izo6mlGjRmGxWMjNzWX58uX069ePadOmMXHiRLZs2UJqairDhw8nISGhro0fV6sUOR/0uCiRX1+UyB0OF7u+OsCnmXb20gHPadeKQ7Y4DhHHWweg694dXBZWQXLf7iT27113A+B85Khx8N036ew48D1f1LSjMKQ9WHvXW4QXIMjjYpi5lHEXxTJ40BAsQW1i1RQRkXOCZqj8pKUzVCdqXPz8nYaVDNu616f0ISK4Zc/nn+0MFUBmZibjx4/HarXy6aefcsEFF7BixQp+//vfM2TIEN577716g7K0tDQ8Hg99+/at186mTZu48847iY+PJyUlpd5/l2nTppGSksLkyZNZvHgxNtsPd9qzs7Nxu9307n2qkpaKUmiGKlC11p3EwvKTbNr5HV8UVJFtbrpcd7jzJANcxQyMgkGJcXRO6oM5tJFn284R7soKDqdlsDe7kK/LzXxrjaGqkSIftS702BnVLYIxw/ucdZGJs+Xvu8zSdik25Ez8HR+aoZLz2pnKp7/yyitcd911dX/v3bs3f/jDH1iwYAGzZ8/m6aef5oknniA8PJznnnuuwR3uiy++uNF2r7zySiZNmsSaNWs4cOBA3XF79uwhJSWF2NjYBskUQK9evbw9TZHzUmy7MKaOHchUIL/kOF/sSmfb0RoyLfUr01VYwthuuYDtVcB+iP1qD5c4CkmKDqJP91i69euDJaphNbtA4Sgp5mBqBpmHikg77uGb4M4cC4kCukNow+NNHjdJrlJGdg1jxPCLiY1KavU+i4icb5RQScA6U9n0xpKtGTNm8Nlnn7FhwwZuuukmTp48yZIlS+o9cne66upqNm3axJ49ezh27Bg1NTXAqdkrgJycnLqEasuWLQDcfPPNDZIpEWmZ+I6R/Oy6YfwMOFJYyhe70th51EG6tRMuc/0Z78LQDqwP7cB6F5ADoZl59KraTZ9QB33iIkns053YC+KxWNvez5+zxsGRg4fJyDxMZuEJMqqDyQmNpSaoPdAeIhv/nNXt4GJnCZd2CeOy5Ivp2KFv4weKiIhPtL1fFGmWcKuZ16f08frz/ipKEW417rl9b8qmL1myhBEjRnD8+HHGjRvH9OnTGz0uLS2NO++8k7y8vCbbOnHiRN2f8/PzAejRo8dZ90lEmq9zbHtunjiSm4HKE5Xs/zaLvXl2vq60kmeJbnB8VVAIqeHdSAUoBoorsWw9QJcaO108FXS1OIm3mYiPCiGuUxRRMR2xdorFZHC5do/HA5UV1BQXUlZYwpGScvLLqsk/6eF7p5UCUzhHgqP/lyDGgCXmjL/QPZ2lDIxwMqhHRy7u14fQsAGG9ldERJpPCVWAMplMLXoXyWI59Vmn+fx6he6TTz6hqqoKOPVeVWVlZYMZJY/Hw6xZs8jLy+OOO+7gjjvuoEePHoSHh2MymVi4cCHPP/98o+/A6X0gkdZji7AxfMQAho849feSE9V8nXaI/XnHyDwBh8yRuE0Nb+I4zRbyQjuRR6dTGxz8L9kCDlQT7kgn2llJlKea6CA30cEQajFjxY0VN8G4sZo8WHEThAcnJpyYqfGYcWDCgZkazJx0erDXeChzBWE3BVNmDafSEsapqab/TTdZaVBA4seCPC4SPMfpHWmmf/eOXJJ0AdE2lZAXEWkrlFDJeSM7O5vHHnsMm83G6NGjWbduHY899hiLFi2qd1xmZiaZmZkMHDiQp59+ukE7hw4darAtPj4egNzcXJ/0XUR+WseIEK4e3oerh5/6+0mHi+ycAjKy88koPkmmI4QjwdE/2U6F1UaF1Ua9sjfeFEQ10WAB3Z/8iMdNvKOMPsHV9IkJp0/vbvTsEUewqvKJiLRZSqjkvOB0OrnvvvuorKxk8eLFTJ48mWuvvZbXX3+dq6++ul4BC7vdDtBoxRe73c7nn3/eYPsVV1zBn//8Z9577z3mzZtHWNiZq4xZraduSbtcrVu2XuR8EmYNot+F3eh3Ybe6bVUVJynILyS/0M73pSfJr3CRX2MmnzCOmxup8uAj7dxVdDWdJD7YTXyEha7tbcTHRdM5PoaQ0Nbrh4iItJwSKjkvLFmyhL179zJ+/Hhuv/12AJ577jkmT57MvHnzGDJkSN0aUD179sRsNrN161ays7PrKvRVVVUxf/78uoTrdIMHD2bkyJFs27aNhx9+mEWLFtVLqnJycnC5XHVl0+Pi4gDIysry5WmLyI+EhofRs08PevZp+L5jjcuNvdKBvbQMe7Edu70c+/GTlFU4qPZQ73E+B0E4/veon6X2McD//bP27yFmE1G2YNq3sxEVHUl0p/ZEd2hHdKgVa5AeDxYROVcooZKANWfOnCb3de3alXnz5gGwc+dO/va3vxETE8PixYvrjhk8eDBz5szhL3/5C3PnzmXlypWYTCY6derEbbfdxuuvv864ceMYNWoUoaGh7Ny5E5fLxdSpU1m9enWD73z22WeZOnUqa9asYePGjSQnJ9ct7JuamsqSJUvqEqoxY8YQGhrKSy+9RHp6OnFxcZhMJmbNmlV3jIi0ruAgM7GRIcRGxkJ3LbItIiLNo4RKAtY777zT5L6+ffsyb948jh8/zuzZs3G5XCxZsoSOHTvWO2727Nls3LiRjRs38uqrr/KrX/0KgIULF5KYmMhbb73F1q1biYyM5IorruB3v/sdb7/9dqPf2aVLF9auXctLL73ERx99xObNm7FYLMTHxzNz5kwuv/zyumM7d+7M8uXLWbZsGTt37qSiogKAW265RQmViIiISAAxeRorVSY+V1BQ0GiVuNN5PB5KSkro2LGj4dXj/FU2XQKDr+PDl7EtvuXvFe2l7VJsSFMUG3Im/o4Pk8nU6HvzZ0Nlg0RERERERLykhEpERERERMRLSqhERERERES8pIRKRERERETES0qoREREREREvKSESkRERERExEtKqERERERERLykhKoNM5lMmEwm3G63v7siYiiXy1UX3yIiIiKBTAlVGxccHExVVZW/uyFiqOrqaoKDg/3dDREREZEWs/i7A3JmYWFhlJWVARASEkJQUJAh7Xo8nnr/FDmdr+LD5XJRXV1NVVUVUVFRhrYtIiIi4g9KqNo4i8VCVFQUJ0+epKyszLABrtl8anJSjxNKY3wVHyaTieDgYKKiorBYdPkRERGRwKcRTQCwWCxERkYCxs0YxMTEAFBUVGRIe3Ju8VV86J0pEREROdcooQowRg1Ia9vRAFcao/gQERERaR4VpRAREREREfGSEioREREREREvKaESERERERHxkhIqERERERERLymhEhERERER8ZKq/PlJW6me1lb6IW2T4kOaotiQpig2pCmKDTkTf8WHEd9r8hi1sJGIiIiIiMh5Ro/8iYiIiIiIeEkJ1Xlq/vz5zJ8/39/dkDZK8SFNUWxIUxQb0hTFhpzJuRAfeofqPFVTU+PvLkgbpviQpig2pCmKDWmKYkPO5FyID81QiYiIiIiIeEkJlYiIiIiIiJeUUImIiIiIiHhJCZWIiIiIiIiXtA6ViIiIiIiIlzRDJSIiIiIi4iUlVCIiIiIiIl5SQiUiIiIiIuIlJVQiIiIiIiJeUkIlIiIiIiLiJSVUIiIiIiIiXlJCJSIiIiIi4iUlVCIiIiIiIl6y+LsDYoyamhref/99tm7dSnFxMREREQwcOJBp06bRsWPHs2qroqKCd955h507d2K324mOjmb48OFMnTqV8PBwH52B+IoRsVFRUcGePXv46quvyM3Npbi4GJPJRLdu3bj88su59tprsVh0OQlERl47TldQUMBDDz2Ew+Fg4MCBLFiwwMBeS2swOjaOHDnC+++/z759+7Db7YSGhtKlSxeSk5O58cYbfXAG4itGxsbevXtZu3YtWVlZVFZWEh4eTu/evZk4cSIDBgzw0RmIL2RnZ/PNN9+QmZlJRkYGpaWlWK1WXn/9da/aC6TxqMnj8Xj83QlpmZqaGp544gnS09Np3749SUlJFBUVkZmZSbt27XjyySfp3Llzs9o6fvw4jzzyCAUFBcTFxdGrVy8OHz5MXl4enTt35qmnniIyMtLHZyRGMSo23nrrLdasWYPJZKJnz5507tyZ8vJy0tPTcTgcJCUlsWDBAkJCQlrhrMQoRl47fuzxxx8nNTUVj8ejhCoAGR0bO3fu5K9//StOp5OEhAS6dOnCiRMnOHToECEhITz33HM+PBsxkpGx8eGHH/LPf/4Tk8nERRddRIcOHTh69ChZWVkAzJw5k2uvvdaXpyMGWrRoEV9++WW9bd4mVIE2HtUt5XPAe++9R3p6OhdeeCGPPPIIoaGhwA8XqhdeeIHHH3+8WW2tWLGCgoICkpOTeeCBBwgKCgJg+fLlrFu3jhUrVnDvvff67FzEWEbFRmhoKDfffDPjx4+nQ4cOddsLCgp44oknOHDgAO+++y633367z85FjGfkteN0n332Gfv372fs2LGsX7/e6G5LKzAyNnJzc3nmmWcICwtj3rx5JCUl1e1zu93k5OT45BzEN4yKjfLyct544w0sFguPPvpovbjYvn07y5YtY+XKlYwePbruO6Rtu/DCC0lISCAxMZHExETuuusur9sKtPGo3qEKcE6nk3Xr1gHw61//ut5FZ9KkSfTo0YO0tDSys7N/si273c6WLVsICgpi5syZdcELcMcdd9CuXTtSUlKw2+2Gn4cYz8jYmDx5Mrfddlu9ZAqgS5cudUnU1q1bDey9+JqR8XG6srIyVq5cyYABAxg1apShfZbWYXRsvPrqqzidTu655556g2YAs9lMYmKicZ0XnzIyNjIyMnA6nfTv379BXIwYMYLu3btTXV3N4cOHjT0J8ZnJkyczdepUhg4dSnR0tNftBOJ4VAlVgDtw4AAVFRXExcXRs2fPBvsvvfRSgAZTsI3Zs2cPHo+Hvn37NvgfwWq1MnToUNxuN3v37jWi6+JjRsbGmSQkJABQWlraonakdfkqPl599VVqamr4zW9+Y0g/pfUZGRuHDx8mLS2NLl26MHToUMP7Kq3LyNiwWq3N+s6IiIiz66QEvEAcjyqhCnAHDx4EaPTCBtCrV696x7Wkrdrtubm5Z9tN8QMjY+NMjh49CtCiu1HS+nwRH7t372bbtm3cfPPNXr97Jf5nZGx8++23AFxyySXU1NSwadMmli9fzvLly9mwYQOVlZUG9Vpag5GxkZiYiM1m49tvv+XAgQP19u3YsYNDhw5x0UUX6VpyHgrE8ajeoQpwxcXFAE1W1al9RKv2uOa09ePHumrVfkdz2hL/MzI2zmTt2rUADBs2rEXtSOsyOj6qqqp45ZVXiI+PZ/LkyYb0UfzDyNjIy8sDIDg4mIcffpj8/Px6+9944w0efPBB+vbt25IuSysxMjbCw8OZNWsWzz77LH/84x/rilIUFhaSlZXFoEGDuOeee4zrvASMQByPKqEKcFVVVQBNVlerfb659riWtFW7vbq6+qz7Ka3PyNhoyieffMK+ffsIDw/XIDrAGB0fb731FkVFRTz66KMqoR/gjIyNiooK4NSNl/DwcB566CH69++P3W7nX//6FykpKSxevJilS5fSvn17g85AfMXo68aIESOIiIhg2bJl9WapoqKi6NevX5uq4iatJxDHo/rVC3A/VfX+bKri1x5rMpla1CdpG4yMjcakpqby2muvYTKZuPvuu5u8kyRtk5HxkZWVxbp16xg9ejT9+/dvadfEz4yMDbfbDYDL5eK+++5j4MCBANhsNmbPnk1BQQFZWVl8/PHHTJ8+3ftOS6sw+nflgw8+YNWqVXVrC8XGxlJYWMjbb7/NqlWryMjI4MEHH2xJlyUABeJ4VO9QBbiwsDCg6Sy9dntzSo7WttXUnaXatrTWUGAwMjZ+7ODBgyxevBin08mMGTNITk72vqPiF0bFh8vl4h//+Ac2m41f/OIXxnZS/MLIa0ftMR06dKhLpk531VVXAbB//36v+iqty8jYSE1NZeXKlSQkJDB37ly6d+9OaGgo3bt358EHH6Rnz57s2LGDr7/+2rgTkIAQiONRzVAFuE6dOgFQUlLS6P5jx47VO645bdV+5sdqv6M5bYn/GRkbpzty5AhPPfUUFRUVTJkyhQkTJrSso+IXRsVHSUkJubm5REdHs3Tp0nr7ah/3yszM5LHHHiM0NJT58+e3tOviY0ZeO2JjYwGIiYlpdH/t9vLy8rPup7Q+I2Nj8+bNwKnKgGZz/fv7ZrOZ5ORkcnJy2L9/f6PJuJy7AnE8qoQqwPXo0QOgyYURa9eCqD2uJW3Vbm9OW+J/RsZGrWPHjvHkk09it9u5/vrrmTJlSss7Kn5hdHzY7fYm1wSpqKggNTUVm8129h2VVmdkbNQuq3DixIlG9x8/fhzwbqZcWp+RsVE7WK6djfix2u1NxY6cuwJxPKqEKsAlJSVhs9k4evQoOTk5DUpM7tixA4AhQ4b8ZFuDBg3CZDKRlpZGWVkZUVFRdfscDgdfffUVJpOJwYMHG3sS4hNGxgac+lF76qmnKCws5Morr+SXv/yl4X2W1mNUfMTGxrJ69epG9+3fv5/HH3+cgQMHsmDBAmM6Lj5n5LVjwIABhISEcOTIEYqLixvcUU5NTQWaLo8sbYuRsVE7xsjKymp0f+322llOOX8E4nhU71AFOIvFwnXXXQfA8uXL6z1v+uGHH3Lw4EGSkpLo3bt33fZ169YxZ84c3njjjXpttW/fnlGjRuF0Onn55ZdxuVx1+1atWkV5eTmXX3651hsKEEbGRnV1NQsXLiQvL4/LLruMWbNmBdTLotKQkfEh5xYjYyMkJIQJEybgcrl4+eWX67W1d+9eNm/ejMlkYuzYsT4+KzGCkbFR++5tSkpKg4WAd+3aRUpKCiaTSe/onsPOpfGoZqjOAbfccgv79u0jPT2d+++/n6SkJIqLi8nIyCAyMrLBOg7l5eXk5+dTWlraoK0ZM2aQkZHBjh07mDNnDomJieTl5ZGXl0dcXJxmJQKMUbHx5ptvkpGRgdlsJigoiBdeeKHR7/vtb3/rs3MR4xl57ZBzi5Gxceutt5KWlsbu3bu5//776d27N+Xl5Xz33Xd4PB6mT59ebwAubZtRsTF8+HBGjBjB9u3bWbRoEYmJicTExFBUVFQ3OzV9+nTi4+Nb7dykZXbv3s27775bb5vT6az3hMLPfvazuhnMc2k8qoTqHBAcHMwf//hH3nvvPVJSUti1axfh4eGMGTOGadOmndVLe+3atWPhwoWsXr2aXbt2sXPnTqKiorjuuuuYOnUqERERPjwTMZpRsVFbXMDtdpOSktLkcUqoAouR1w45txgZG7Vt/ec//2HLli3s3bsXq9VK//79mThxYrMfO5a2wajYMJlMPPDAA2zcuJHNmzdz6NAhcnNzsdlsDB48mAkTJjBo0CDfnowYqry8nIyMjHrbPB5PvW3NLUATaONRk6eli9GIiIiIiIicp/QOlYiIiIiIiJeUUImIiIiIiHhJCZWIiIiIiIiXlFCJiIiIiIh4SQmViIiIiIiIl5RQiYiIiIiIeEkJlYiIiIiIiJeUUImIiIiIiHhJCZWIiIiIiIiXlFCJiIiIiIh4SQmViIiIiIiIl5RQiYiIiIiIeEkJlYiIiIiIiJeUUImIiIiIiHhJCZWIiIiIiIiXlFCJiIiIiIh4SQmViIiIiIiIl/4/2jz6jxeVTb8AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 960x720 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "plt.style.use('ggplot')\n",
    "plt.rcParams['figure.dpi'] = 150\n",
    "\n",
    "plt.plot(x[1:-1], u)\n",
    "\n",
    "xe = np.linspace(0, 1, 101)\n",
    "plt.plot(xe, -0.5*xe + 0.5*xe**2)\n",
    "plt.legend([\"Computed\", \"Exact\"])"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.13.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}